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Abstract 

The  US  Air  Force  uses  many  combat  simulation  models  to  assist  them  in 
performing  combat  analyses.  BRAWLER  is  a  high-resolution  air-to-air  combat 
simulation  model  used  for  engagement-level  analyses  of  few-on-few  air  combat. 
THUNDER  is  a  low-resolution  combat  simulation  model  used  for  campaign-level 
analyses  of  theater-level  warfare.  BRAWLER  is  frequently  used  to  ensure  that 
THUNDER  air-to-air  inputs  are  valid.  This  thesis  describes  the  confederation  of 
THUNDER  and  BRAWLER  by  clearly  showing  how  one  particular  BRAWLER  output, 
the  effectiveness  of  a  missile  type,  is  transformed  into  THUNDER  air-to-air  input  data 
Since  BRAWLER  is  a  stochastic  simulation  model,  it  is  necessary  to  replicate  a  number 
of  BRAWLER  simulation  runs  in  order  to  obtain  a  sufficiently  accurate  estimate  of  the 
mean  missile  effectiveness,  a  number  that  varies  for  each  different  BRAWLER  combat 
scenario.  This  thesis  focuses  on  using  two  different  sequential  methods  to  determine 
when  the  minimum  number  of  BRAWLER  runs  has  been  performed  to  obtain  a  specified 
relative  precision.  One  method  uses  classical  statistical  analysis  techniques,  while  the 
other  uses  the  more  modern  technique  of  bootstrap  resampling.  The  performance  of  these 
two  methods  is  compared. 
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CROSS-RESOLUTION  COMBAT  MODEL  CALIBRATION 
USING  BOOTSTRAP  SAMPLING 

1.  Introduction 


1.1.  The  Simulation  and  Analysis  Facility 

The  Simulation  and  Analysis  Facility  (SIMAF)  is  a  $7  million  center  currently 
being  developed  by  the  United  States  Air  Force’s  (USAF)  Aeronautical  Systems  Center 
(ASC).  The  SIMAF  vision  is  to  become  the  premiere  facility  for  all  ASC  modeling, 
simulation,  and  analysis  (MS&A)  to  support  decisions  concerning  the  acquisition  of 
weapon  systems  at  ASC.  The  SIMAF  program  is  managed  by  ASC/SM.  In  May  1997, 
ASC/SM  requested  that  the  Air  Force  Institute  of  Technology  (AFTT)  provide  them  with 
technical  guidance  and  analysis  on  a  pilot  effort  to  demonstrate  SIMAF's  capability.  In 
particular,  SIMAF  wanted  AFTT  to  help  produce  a  proof-of-concept  scenario 
demonstrating  the  utility  of  simulation  environments  in  the  acquisition  cycle. 

The  SIMAF  operational  concept  is  that  the  center  will  be  a  "virtual  facility"  as 
well  as  a  "physical  facility."  As  a  virtual  facility,  it  will  function  as  an  integrator  of 
products,  tools,  and  technology  and  as  a  communications  and  networking  node.  As  a 
physical  facility,  it  will  function  as  a  physical  gateway  for  MS&A  to  a  synthetic 
battlespace  environment,  provide  a  common  environment  for  both  constructive  and 
virtual  analysis,  and  support  multiple  weapon  acquisition  programs  at  different  levels  of 
security  (ASC/SM,  1997). 
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1 .2.  SIMAF  Proof  of  Concept 


To  demonstrate  proof  of  the  SIMAF  concept,  a  pilot  program  was  created  to 
demonstrate  the  so-called  cross-resolution  modeling  process;  specifically,  the 
confederation  of  an  engagement-level  simulation  and  a  campaign-level  simulation.  The 
campaign-level  model  used  in  this  program  is  THUNDER,  and  the  engagement-level 
model  is  BRAWLER.  THUNDER  is  the  USAF’s  primary  campaign  analysis  tool,  and 
BRAWLER  is  its  primary  air-to-air  engagement  analysis  tool.  Both  THUNDER  and 
BRAWLER  are  considered  legacy  models.  BRAWLER  became  operational  in  1976, 
while  THUNDER’S  first  production  release  was  ten  years  later  in  1986.  While  the 
current  trend  in  modeling  and  simulation  is  toward  virtual  models,  both  of  these  legacy 
models  are  constructive  models.  These  two  particular  models  were  chosen  for  three 
reasons: 

1 .  There  is  a  wealth  of  information  available  about  each. 

2.  There  are  many  expert  users  of  the  models  available  for  consultation. 

3.  Both  models  are  part  of  the  USAF  Analysis  Toolkit. 

As  part  of  the  USAF  Analysis  Toolkit,  both  will  continue  in  use  while  the  USAF 
transitions  from  current  legacy  models  to  new  full-dimensional  models  (Langbehn,  1998). 

An  effort  was  made  to  include  a  virtual  engagement-level  model  in  this  research, 
namely  the  USAF’s  second-generation  version  of  the  Man-in-the-Loop  Advanced  Air-to- 
Air  System  Performance  Evaluation  Model  (MIL-AASPEM II).  However,  the  newness 
of  this  model  precluded  its  use,  since  at  the  beginning  of  this  research  undertaking  neither 
a  production  version  of  the  model  nor  appropriate  training  was  available. 
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1.3.  Cross-Resolution  Modeling 


This  research  demonstrates  the  confederation  of  THUNDER  and  BRAWLER. 
BRAWLER  is  a  high-resolution  model  and  THUNDER  is  a  low-resolution  model. 
THUNDER  models  many  other  aspects  of  warfare  besides  air  combat,  so  in  order  to 
avoid  being  a  prohibitively  large  model,  its  air-to-air  model  is  less  detailed  than 
BRAWLER’S.  Their  confederation  is  important  because  BRAWLER  is  used  to  ensure 
that  THUNDER  air-to-air  inputs  are  valid,  or  in  other  words,  to  calibrate  THUNDER. 

We  say  that  THUNDER  input  data  that  has  been  properly  calibrated  by  BRAWLER  has  a 
pedigree. 

This  research  clearly  shows  how  a  particular  BRAWLER  output  metric,  called  a 
measure  of  effectiveness  (MOE),  is  transformed  into  input  data  for  THUNDER.  The 
specific  BRAWLER  MOE  used  in  this  research  is  a  measure  of  missile  effectiveness 
called  adjusted  kills  per  firing  (AKPF).  The  AKPF  MOE  is  calculated  from  BRAWLER 
missile  effectiveness  data  using  a  process  introduced  in  Section  3.4.5. 

When  we  calibrate  THUNDER  using  BRAWLER,  we  might  say  that  this  missile 
effectiveness  data  crosses  a  level  of  model  resolution  when  passing  from  BRAWLER  to 
THUNDER.  By  using  this  calibration  process,  we  follow  an  important  principle  of  cross¬ 
resolution  modeling  is,  which  is  to  ensure  that  a  pedigree  for  the  data  is  established  before 
crossing  a  level  of  resolution.  In  this  calibration  process,  THUNDER  input  data  has  a 
direct  lineage  to  BRAWLER  output  data,  and  is  therefore  called  a  pedigree  database. 
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1 .4.  Stochastic  Simulation 


In  the  current  THUNDER  calibration  process,  the  AKPF  MOE  from  BRAWLER 
is  used  to  calibrate  THUNDER  input  for  air-to-air  probability  of  kill  (PK),  which  deals 
with  the  likelihood  that  an  aircraft  with  a  given  weapon  will  kill  its  target.  However,  the 
AKPF  may  vary  drastically  for  each  realization  of  a  BRAWLER  engagement.  This  is 
because  BRAWLER  is  a  stochastic  simulation  model  involving  many  random  variables. 
This  requires  that  BRAWLER  engagements  be  replicated  multiple  times.  The  mean  of 
the  AKPF  values  obtained  from  the  replications  of  an  engagement  is  our  best  estimate  of 
the  expected  value  of  AKPF  over  the  long  run,  which  we  call  the  expected  AKPF.  To 
determine  if  our  observed  mean  AKPF  is  sufficiently  accurate  for  this  purpose,  a  standard 
confidence  interval  is  often  calculated.  If  a  sufficiently  narrow  confidence  interval  about 
the  mean  is  obtained,  we  are  confident  in  our  estimate  of  the  true  AKPF  is  a  good  one. 
Therefore,  we  are  confident  of  providing  THUNDER  with  a  good  PK  input  value. 

The  number  of  replications  required  to  obtain  a  sufficiently  accurate  estimate  of 
the  true  AKPF  varies  for  each  BRAWLER  scenario.  A  BRAWLER  scenario  is  defined 
by  information  such  as  the  types  and  numbers  of  aircraft  involved  in  combat  and  the 
length  of  the  simulation  run.  Different  combinations  of  the  numbers  of  aircraft  involved 
in  combat  are  referred  to  as  engagement  force  structures.  An  example  of  an  engagement 
force  structure  is  four  friendly  aircraft  versus  four  hostile  aircraft  (4v4).  To  obtain  an 
accurate  estimate  of  the  true  AKPF,  we  must  ensure  that  this  mean  is  calculated  from  a 
number  of  BRAWLER  replications  that  is  greater  than  or  equal  to  N,  where  N  is  the 
number  of  replications  necessary  to  meet  a  precision  criterion.  One  such  criterion  that  is 


1-4 


often  used  by  BRAWLER  users  to  estimate  means  of  MOE’s  and  used  in  this  thesis  is  the 
follows: 

•  Enough  replications  have  been  performed  when  the  endpoints  of  a 

100  •  (1  -a)%  standard  confidence  interval  about  the  estimated  mean  are 
within  ±y%  of  the  estimated  mean. 

Here  a  is  the  significance  level.  Another  way  of  stating  this  criterion  is  that  the  interval 
satisfies  the  requirement  that  the  relative  error  e  of  the  estimated  mean  is  equal  to  y . 

This  interval  about  the  estimated  mean  is  an  approximate  100-(l-a)%  confidence 

interval  for  the  expected  AKPF.  We  may  also  consider  this  precision  criterion  as  a 
stopping  criterion  for  performing  BRAWLER  replications. 

1.5.  The  Campaign  Analysis  Group 

One  USAF  analysis  group  that  uses  BRAWLER  and  THUNDER  is  the 
Aeronautical  Systems  Center  Development  Planning  Campaign  Analysis  Group 
(ASC/XRA),  who  are  a  major  potential  user  of  SIMAF.  ASC/XRA  performs  campaign- 
level  analyses  for  the  USAF’s  top-priority  weapon  system  programs  such  as  the  Joint 
Strike  Fighter  (JSF)  program.  The  JSF  is  currently  being  developed  to  replace  the  F-16 
and  perhaps  other  aircraft  as  well.  Over  the  past  several  years,  ASC/XRA  has  performed 
many  BRAWLER  replications  of  many  different  scenarios  and  has  calculated 
100  •  (1  -  a)%  standard  confidence  intervals  for  many  different  combinations  of  the 
following: 

1 .  the  number  of  BRAWLER  replications  performed, 

2.  the  value  of  the  significance  level  a ,  and 

3.  the  type  of  BRAWLER  scenario  that  is  performed. 
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These  confidence  intervals  were  used  to  help  estimate  how  many  BRAWLER  replications 
of  a  scenario  are  actually  necessary  to  meet  their  precision  criterion,  and  these  values  are 
currently  used  as  guidelines  for  similar  scenarios.  BRAWLER  users  in  ASC/XRA  use 
these  guidelines  as  well  as  their  own  experience  to  choose  the  number  of  replications  they 
wish  to  perform.  For  many  4v4  engagements,  ASC/XRA  typically  performs  between  200 
and  300  BRAWLER  replications  (Taranto,  1998). 

This  information  brought  the  following  questions  to  mind: 

1 .  Could  another  approximate  (1  -  a)%  confidence  interval  besides  the  standard 

approximate  confidence  interval  be  calculated  such  that  the  stopping  criterion 
is  met  with  fewer  replications? 

2.  If  so,  how  would  this  interval  compare  with  the  standard  confidence  interval 
in  its  coverage  of  the  expected  AKPF? 

1.6.  Bootstrap  Sampling 

It  turns  out  that  one  such  method  for  calculating  such  an  interval  involves  using  a 
technique  called  bootstrap  sampling,  introduced  in  Section  4.5.  The  bootstrap  sampling 
method  involves  analyzing  pseudo-data  obtained  by  resampling,  or  “bootstrapping”,  the 
original  data  a  large  number  B  of  times.  Bootstrapping  is  a  technique  recommended  by 
Efron  (1982)  for  constructing  approximate  confidence  intervals.  Using  bootstrap 
sampling,  we  can  calculate  a  bootstrap  percentile  interval.  A  (1  -a)%  bootstrap 
percentile  interval  is  another  approximate  (1  -a)%  approximate  confidence  interval  for 
the  expected  AKPF,  and  has  subtle  differences  from  a  standard  confidence  interval. 
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1.7.  Purpose  of  Research 


The  main  purpose  of  this  research  is  to  investigate  if  the  number  of  BRAWLER 
replications  necessary  to  obtain  a  (l-a)%  bootstrap  percentile  interval  that  meets  the 
stopping  criterion  above  is  less  than  the  number  required  to  obtain  a  (1  -a)%  standard 
confidence  interval  that  meets  the  same  criterion,  where  in  both  cases  a  =  y  =  .10 .  In 

this  investigation,  200  BRAWLER  replications  were  performed  for  each  of  three  different 
scenarios,  which  are  described  in  Section  4.2.4.  The  specific  objective  of  the 
investigation  can  be  described  in  four  steps  to  complete  for  each  scenario,  as  follows: 

1.  Calculate  the  estimated  mean  AKPF  and  a  standard  90%  confidence  interval 
about  the  estimated  mean  using  all  200  BRAWLER  replications.  Use  this 
mean  and  confidence  interval  as  our  best  estimate  of  the  expected  value  of 
AKPF.  We  call  this  estimate  our  Truth  model ,  since  we  are  using  it  to 
represent  the  true  expected  AKPF. 

2.  Determine  the  minimum  number  IV,  of  BRAWLER  replications  necessary  to 
provide  an  approximate  90%  standard  confidence  interval  such  that  the 
relative  error  of  the  estimated  mean  is  equal  to  .10.  Use  this  estimated  mean 
as  another  estimate  of  the  expected  value  of  AKPF.  A  sequential  procedure 
that  accomplishes  this,  denoted  Candidate  SI,  is  described  in  Chapter  4.6.2. 

3.  By  bootstrapping  the  AKPF  data  B  times,  where  B  =  501 ,  determine  the 
minimum  number  AT2  of  BRAWLER  replications  of  each  scenario  necessary 
to  provide  a  bootstrap  90%  percentile  interval  that  meets  the  same  criterion  as 
in  step  2.  Use  this  estimated  mean  as  yet  another  estimate  of  the  expected 
value  of  AKPF.  A  sequential  procedure  that  accomplishes  this,  denoted 
Candidate  S2,  is  described  in  Chapter  4.6.4. 

4.  Compare  AT,  and  AT2  and  the  two  intervals  generated  in  steps  2  and  3. 

Answer  the  following  questions: 

a.  Is  the  expected  value  of  N2  less  than  the  expected  value  of  TV,  (denoted 
E(TV2)  <  E(  AT, ))?  In  other  words,  can  90%  bootstrap  percentile  intervals 
using  N2  replications,  where  AT2<  AT, ,  be  consistently  constructed  that 
meet  the  same  stopping  criterion  as  a  standard  90%  confidence  interval 
using  AT,  replications? 
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b.  Do  either  (or  both)  of  the  two  intervals  calculated  in  steps  2  and  3  overlap 
the  interval  calculated  in  step  1?  In  other  words,  are  the  corresponding 
estimates  of  AKPF  statistically  similar? 

The  procedures  used  in  Candidates  SI  and  S2  are  displayed  in  Figure  1.1. 


Candidate  1 : 

After  each  run 
following  n0  pilot  runs 


n0  pilot  runs 


..•gat*! 


at  N1  when 
Relative  Error 
£  <  .10 


Candidate  2: 

After  each  run  following 
nQ  pilot  runs, 


Bootstrap  current  n  data 
points  B=501  times,  then... 


n0  pilot  runs 


(□□□□□□ 


□  □□□□□  \  bootstrap 
°  ?  °  °  °  °  /  B=  501 


Calculate  90% 
Percentile  Interval  (PI)  for 
B  =  501  Bootstrapped  Means 

Median  of 

Bootstrapped  Means 


at  N2  when 
Relative  Error 
e<.10 


!□□□□□□/ 


Figure  1.1.  Procedures  for  Candidates  1  and  2 


1.8.  Thesis  Format 

Chapter  2  contains  the  results  of  a  review  of  the  literature  on  several  different 
subjects  related  to  this  thesis.  Chapter  3  provides  a  background  discussion  on  the  areas 
that  are  explored  in  this  thesis.  Chapter  4  discusses  the  methodology  used  to  obtain  the 
results  of  this  research,  and  Chapter  5  presents  those  results.  Chapter  6  summarizes  the 
research,  discusses  the  results,  and  presents  ideas  for  additional  research  in  related  areas. 
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Supporting  material  for  this  thesis  are  contained  in  the  appendices.  Appendix  A 
contains  a  Mathcad  6.0  program  written  by  the  author  and  designed  to  run  on  a 
Windows™  personal  computer  (PC)  called  bootstrap.mcd  that  performs  bootstrapping 
functions  and  contains  algorithms  that  determine  the  number  of  BRAWLER  replications 
required  to  accurately  estimate  mean  AKPF.  Appendix  B  contains  a  FORTRAN  77 
program  written  by  Johnson  and  Lucas  (1994)  called  enoughruns.F  that  is  designed  to 
run  between  BRAWLER  replications  on  a  Sun  or  Silicon  Graphics  workstation.  It  uses 
bootstrapping  to  determine  if  enough  BRAWLER  replications  have  been  performed  to 
accurately  estimate  the  expected  value  of  an  MOE  called  Exchange  Ratio.  Appendix  C 
contains  a  BRAWLER  post-processing  A  WAT  program  called  STATS  written  by  Taranto 
(1997).  STATS  calculates  many  BRAWLER  MOE’s  using  BRAWLER  replication  log 
files  and  generates  an  output  file  called  STATS.prt  that  summarizes  these  MOE’s. 
Example  missile  effectiveness  reports  from  STATS.prt  appear  in  Appendix  D. 

Appendix  E  contains  an  AWK  program  called  msl_summary.awk,  a  subprogram  of 
STATS,  in  which  the  AKPF  MOE  is  calculated.  Appendix  F  contains  a  Unix  C-Shell 
program  written  by  the  author  called  make  JlyAKPF  that  runs  STATS  on  each  individual 
BRAWLER  replication  log  file  and  creates  one  fde  of  AKPF  data  that  may  be  used  for 
bootstrap  sampling.  Appendix  G  contains  the  BRAWLER  AKPF  and  Firings  data  from 
200  replications  of  each  of  five  different  BRAWLER  scenarios.  Finally,  Appendix  H 
presents  the  derivation  of  an  equation  for  AKPF. 
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2.  Literature  Review 


2.1.  BRAWLER 

Johnson  (1994)  investigated  how  many  BRAWLER  replications  were  required  to 
obtain  a  relative  error  of  .10  for  a  BRAWLER  measure  of  outcome  (MOO)  called 
exchange  ratio.  Exchange  ratio  is  defined  as  the  number  of  Red  aircraft  killed  (RK)  in  an 
engagement  divided  by  the  number  of  Blue  aircraft  killed  (BK).  The  exchange  ratio  for  a 
single  BRAWLER  replication  i  is  expressed  as  follows: 


ER.. 


RK , 
BK: 


(2.1) 


However,  since  BRAWLER  is  a  stochastic  simulation  model  (see  Section  3.2.3  for 
details),  we  are  more  interested  in  the  mean  exchange  ratio  (ER)  across  A  BRAWLER 
replications.  ER  is  calculated  as  the  sum  of  all  Red  aircraft  killed  ( RK)  divided  by  the 
sum  of  Blue  aircraft  killed  (BK),  defined  as  follows: 


N 

=  f - 

(2.2) 

i= 1 

where  N  is  the  total  number  of  BRAWLER  replications  that  have  been  performed. 

Johnson  found  that  the  number  of  replications  needed  varies  dramatically, 
depending  on  the  number  of  aircraft  in  the  engagement  and  the  resulting  ER.  For 
example,  a  4v4  engagement  with  a  Blue  PK  =  1.00  and  Red  PK  =  0.25,  1 17  runs  were 
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required,  resulting  in  ER  =  2.49.  A  similar  2v2  engagement  required  approximately  400 
runs  resulting  in  an  approximate  ER  of  8.5. 

Taranto  (1995)  describes  a  process  for  calibrating  THUNDER  with  BRAWLER 
developed  by  ASC/XRA.  This  process,  along  with  a  pedigree  database,  as  defined  in 
Section  1.3,  was  developed  for  the  Department  of  Defense’s  (DOD’s)  Joint  Advanced 
Strike  Technology  (JAST)  program,  which  evolved  into  the  JSF  program.  In  this 
calibration  process,  the  performance  of  an  individual  aircraft/weapon  combination  versus 
an  opponent  aircraft  is  used  to  calibrate  THUNDER  air-to-air  engagements  with 
BRAWLER  engagements. 

2.2.  THUNDER 

Webb  (1995)  examined  THUNDER'S  variability  using  a  univariate  analysis.  He 
also  examined  the  model’s  interrelationships  using  a  multivariate  analysis  (principal 
components  analysis  and  factor  analysis).  Finally,  he  examined  THUNDER’S  sensitivity 
to  input  parameters. 

Farmer  (1996)  modified  the  appropriate  THUNDER  input  variables  and  used 
Design  of  Experiments  (DOE)  and  Response  Surface  Methodology  (RSM)  to  quickly 
show  the  effects  of  changing  force  structure.  He  also  used  factor  analysis  to  find 
relationships  among  different  MOE’s  to  create  new,  simpler  variables.  The  methods  he 
used  provide  accurate,  "quick  turn"  analysis  tools  for  a  force  structure  decision-maker  to 
use. 

ASC/XRA  (1996)  ran  several  experiments  using  THUNDER  and  identified  four 
primary  inputs  that  determine  the  outcome  of  THUNDER  air-to-air  engagements.  They 
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determined  that  it  is  possible  to  correlate  two  of  these  four  inputs  with  BRAWLER 
results.  This  resulted  in  development  of  the  methodology  for  THUNDER/BRAWLER 
calibration  mentioned  in  Section  2.1.  In  the  current  calibration  process,  the  BRAWLER 
MOE  of  adjusted  kills  per  firing  (AKPF)  for  a  given  type  of  missile  in  an  engagement 
becomes  the  THUNDER  air-to-air  probability  of  kill  (PK)  input  value. 

Grier  (1997)  developed  another  “quick  turn”  tool  designed  to  identify  the  cost  and 
capabilities  of  alternative  force  structures  using  factor  analysis  and  response  surface 
methodology.  Using  these  techniques,  he  was  able  to  successfully  “link  procurement 
dollars  to  an  alternative  force  structure’s  combat  capability.”  His  work  allows  the  USAF 
to  evaluate,  within  48  hours,  various  force  structures’  combat  capability  in  terms  of 
theater  level  campaign  objectives. 

2.3.  Cross-Resolution  Modeling 

Davis  (1995)  describes  some  of  the  challenges  in  connecting  models  that  were 
developed  independently.  He  recommends  that  connection  activities  for  these  types  of 
models  be  guided  by  design  work  to  identify  how  the  models  would  have  been  developed 
if  designed  with  integration  in  mind.  He  claims  that  this  approach  makes  it  possible  to 
connect  the  models  usefully  and  comprehensively. 

Hillestad,  Owen,  and  Blumenthal  (1995)  examined  differences  in  outcomes 
predicted  by  different-resolution  models  using  identical  combat  situations,  which 
provided  insights  into  the  problems  of  aggregation.  They  found  that  intuition  regarding 
outcomes,  causes,  and  effects  were  frequently  wrong,  leading  to  bad  aggregate 
approximations.  They  note  that  although  scaling  for  different  levels  of  resolution  is 
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possible,  they  were  unable  to  find  a  method  of  predicting  the  appropriate  resolution  level 
scaling  technique  and  factors  in  advance. 

Hillestad  and  Moore  (1996)  investigated  and  implemented  alternatives  to  different 
modeling  aspects  they  believe  are  essential  to  improving  theater-level  campaign  analysis. 
These  aspects  include  how  to  create  more  flexible  structures  to  simulate  the  wide  range  of 
future  scenarios  and  their  uncertainties,  and  how  to  link  low-resolution  models  to  high- 
resolution  models  in  an  analytically  valid  way.  They  built  the  Theater-Level  Campaign 
(TLC)  model  to  improve  their  capability  to  perform  analysis  in  the  future. 

2.4.  Bootstrap  Sampling  Method 

Efron  (1979)  is  the  author  of  the  first  largely  read  introduction  to  the  bootstrap 
and  is  often  credited  with  popularizing  the  bootstrap  as  a  way  to  estimate  a  sampling 
distribution.  Efron  (1982)  presents  more  details  about  the  bootstrap  and  its  relationship 
to  more  well  known  methods  such  as  the  jackknife.  Efron  and  Tibshirani  (1993)  offer  a 
comprehensive  review  of  bootstrap  resampling  techniques. 

Johnson  and  Lucas  (1995)  developed  a  FORTRAN  77  program  called  enoughruns, 
appearing  in  Appendix  B,  that  uses  the  bootstrap  percentile  interval  method,  introduced  in 
Section  4.5. 3.4.  The  program  enoughruns  is  focused  specifically  on  the  exchange  ratio 
MOO  and  can  be  automatically  called  after  completing  20  BRAWLER  replications,  or 
“runs.”  It  determines  if  enough  BRAWLER  runs  have  been  performed  to  satisfy  a 
specific  stopping  criterion  discussed  in  Section  4.6.3.  If  so,  BRAWLER  replications 
cease;  if  not,  another  BRAWLER  replication  is  performed  and  the  process  is  repeated 
after  the  next  BRAWLER  replication  is  complete.  The  stopping  criterion  is  based  upon 
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the  use  of  a  bootstrap  percentile  interval,  introduced  in  Section  4.5.3.3.  Johnson  and 
Lucas’  methodology  is  discussed  in  detail  in  Section  4.6.3 

Kim  (1996)  used  bootstrap  techniques  in  his  work  on  the  C-17  aircraft’s  personnel 
airdrop  problem.  He  developed  cumulative  distribution  functions  of  the  maximum 
parachute  centerline  entanglement  risk  for  the  C-17  using  bootstrap  techniques,  and  then 
compared  the  effects  of  different  C-17  aircraft  configurations  on  the  entanglement 
distribution  function.  Using  these  techniques,  he  showed  that  certain  configurations  of 
the  C-17  showed  a  lesser  entanglement  risk  than  did  the  older  C-141  aircraft. 

DiCiccio  and  Efron  (1996)  survey  bootstrap  methods  for  producing  approximate 
confidence  intervals,  with  the  goal  of  improving  upon  the  accuracy  of  standard 
confidence  intervals  by  an  order  of  magnitude.  They  also  provide  an  overview  of  four 
different  bootstrap  confidence  interval  procedures. 
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3.  Background 


3.1.  Introduction 

In  this  section  we  provide  a  background  of  topics  pertinent  to  this  research.  We  start  with 
basic  modeling  and  simulation  concepts  in  Section  3.2.  In  Section  3.3  we  discuss  cross¬ 
resolution  modeling  and  some  of  the  issues  involved  with  it.  Finally,  in  Section  3.4  we 
discuss  the  combat  models  used  in  this  research,  BRAWLER  and  THUNDER. 

3.2.  Modeling  and  Simulation 

3.2.1.  Modeling  and  Simulation  Definitions 

The  DOD  defines  constructive  models  as  “a  class  of  models  or  simulations 
executed  by  computer  software  without  man-in-the-loop,”  and  virtual  models  or 
simulations  as  “a  class  of  models  or  simulations  that  feature  live  personnel  (man-in-the- 
loop)  in  a  computerized  simulation  environment.”  It  states  that  the  analysis  of 
engagement-level  models  is  “concerned  with  estimating  the  effectiveness  of  systems 
within  various  classes  of  engagements:  air-to-air,  surface-to-air,  and  air-to-surface  with 
stated  levels  of  performance,”  and  that  “the  resulting  estimates  are  called  Measures  of 
Effectiveness  (MOE’s) .”  It  further  states  that  the  analysis  of  campaign-level  models  is 
“concerned  with  the  cumulative  long-term  effects  of  kills  and  losses  on  the  outcome  of 
theater  level  conflict  of  campaign  duration,”  and  that  it  is  “used  to  determine  how  forces 
should  be  used  to  achieve  campaign  objectives.”  It  adds  that  “outcomes,  measured  in 
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units  that  bear  some  relationship  to  the  overall  objectives  specified  in  the  planning 
scenarios,  are  called  Measures  of  Outcome  (MOO’s)”  (ASC/XRB,  1996). 

Constructive  simulation  models  are  commonly  used  for  engagement-level 
simulations  and  especially  for  theater-level  simulations.  Virtual  models  are  not  currently 
used  at  the  theater  level,  due  to  their  complexity,  but  are  used  effectively  at  lower  levels 
such  as  the  engagement  level.  Both  constructive  and  virtual  simulation  models  are 
common  for  engagement-level  simulations,  with  the  current  trend  towards  virtual 
simulation. 

3.2.2.  Simulation  Model  Resolution 

Campaign-level  simulation  models  are  not  normally  operated  at  the  entity  level; 
rather,  they  are  usually  operated  using  aggregate  groups  of  entities  to  represent  a  single 
object.  They  are  usually  self-contained  and  may  be  very  detailed  in  their  particular  areas 
of  interest,  but  tend  to  be  nominal  in  most  areas  (Hardy  and  Healy,  1994).  For  this 
reason,  campaign-level  models  are  usually  considered  low-resolution  models. 

Engagement-level  simulations  are  usually  modeled  based  on  individual  entities. 
Because  they  are  modeled  this  way,  they  are  usually  far  more  detailed  than  campaign- 
level  simulations  and  are  usually  considered  high-resolution  models.  Because  of 
differences  in  model  levels  and  resolutions,  it  is  useful  to  use  a  set  of  models  to  analyze 
complex  processes  such  as  warfare.  The  core  model  set  that  ASC/XRA  uses  and  the 
hierarchy  of  these  models  is  displayed  in  Figure  3.1.  Recall  that  the  models  we  are 
focusing  on  in  this  research  are  BRAWLER  and  THUNDER. 
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The  main  difference  between  the  methodologies  of  THUNDER  and  BRAWLER 
is  the  level  of  detail  in  which  they  model  air  warfare.  Whereas  air-to-air  attrition  is  very 
detailed  in  BRAWLER,  in  THUNDER  it  is  modeled  at  a  very  aggregated  level.  For 
example,  in  THUNDER,  the  average  effectiveness  of  a  missile  system  versus  a  target 
aircraft  is  approximately  the  same  regardless  of  whether  or  not  the  missile  is  the  primary 
munition,  is  typically  employed  with  specific  tactics,  or  encounters  the  same  threat 
system  with  different  weapons.  In  contrast,  weapons  for  both  Red  and  Blue  are  explicitly 
modeled  in  BRAWLER  and  their  performance  is  verified  against  flight  tests  and 
laboratory  experiments.  Tactics  are  modeled  explicitly  and  are  based  upon  intelligence 
estimates  of  threat  weapon  system  capabilities,  threat  cultural  biases  and  pilot  proficiency 
(ASC/XRA,  1996). 

3.2.3.  Stochastic  Simulation  Models 

Both  BRAWLER  and  THUNDER  are  stochastic  simulation  models  since  they 
involve  many  random  variables.  Each  realization  of  a  BRAWLER  engagement  or  a 
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THUNDER  campaign  is  called  a.  Each  BRAWLER  and  THUNDER  output  must  be 
considered  a  random  variable,  the  value  of  which  varies  from  one  replication  to  the  next 
and  follows  a  specific  statistical  distribution.  Therefore  a  statistical  analysis  of  the  output 
data  is  required,  requiring  many  BRAWLER  and  THUNDER  replications.  As  mentioned 
in  Section  1 .4,  the  AKPF  may  vary  drastically  for  each  replication  of  a  BRAWLER 
engagement. 

A  question  with  stochastic  models  is  “how  many  replications  are  enough?”  When 
trying  to  obtain  a  specified  precision  for  a  statistic  such  as  the  mean  AKPF,  the  answer  to 
this  question  is  not  obvious.  It  varies  for  each  simulation  model,  and  even  for  different 
statistics  of  interest  within  the  same  model.  Law  and  Kelton  (1991)  present  a  sequential 
procedure  for  obtaining  an  estimate  of  the  mean  with  a  specified  relative  error  and 
confidence  level  that  allows  the  simulation  user  to  only  perform  as  many  replications  as 
necessary.  This  procedure  was  used  in  this  research  and  is  discussed  in  Section  4.6.1. 
Johnson  and  Lucas  (1994)  present  a  different  sequential  procedure  that  uses  the  bootstrap 
method,  also  used  in  this  research  and  discussed  in  Section  4.6.3.  The  latter  procedure  is 
contained  in  the  program  enoughruns,  previously  mentioned  in  Sections  2.1  and  2.4. 

3.3.  Cross-Resolution  Modeling 

In  campaign-level  models  such  as  THUNDER,  much  of  the  force-effectiveness 
input  data  such  as  attrition  data  comes  from  higher  resolution  models.  The  linking  of 
high-resolution  models  to  low-resolution  models  is  called  cross-resolution  modeling. 
When  specific  models  are  linked  in  this  manner,  they  are  cross-coupled.  When  models  of 
different  levels  are  linked  together,  they  are  vertically  integrated.  Cross-resolution 
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modeling  and  vertical  integration  is  essential  for  performing  quality  campaign-level 
analyses;  Hillestad  and  Moore  (1996)  list  several  reasons  why,  including: 

1 .  Weapons  system  test  results  are  usually  provided  with  too  much  resolution  for 
a  campaign  model  and  may  cover  only  a  small  subset  of  the  cases  encountered 
in  a  theater  campaign 

2.  Some  phenomena  must  be  modeled  at  higher  levels  of  resolution  to  be 
represented  judiciously  in  theater  models. 

3.  It  is  difficult  to  compare  the  effects  of  weapon  systems  on  campaign  outcomes 
unless  higher-resolution  data  about  those  systems  are  available. 

4.  Detailed  models  are  needed  to  evaluate  new  systems  whose  technological 
characteristics  can  be  defined  but  for  which  test  data  are  not  available. 

Before  beginning  the  cross-coupling  process,  one  must  first  identify  the  high- 
resolution  model(s)  to  use  to  provide  attrition  data  to  the  low-resolution  campaign  model. 
BRAWLER  is  the  air-to-air  combat  model  that  the  USAF  uses  to  provide  air-to-air 
combat  data  to  THUNDER.  In  other  words,  THUNDER  is  cross-coupled  with 
BRAWLER  with  respect  to  air-to-air  engagements.  We  may  also  say  that  THUNDER 
and  BRAWLER  are  vertically  integrated.  We  discuss  in  Section  3.4.5  how  we  calibrate 
the  air-to-air  aspects  of  the  higher-level  THUNDER  model  with  the  lower-level 
BRAWLER  model  using  cross-resolution  modeling  techniques.  BRAWLER  may  also  be 
calibrated  with  still-lower-level  engineering  models  that  perform  specific  aircraft  and 
weapon  functions  in  much  greater  detail  than  BRAWLER  does. 

The  cross-coupling  process  begins  by  generating  a  limited  set  of  inputs  over  a 
range  of  situations  that  are  likely  to  be  encountered  in  the  campaign-level  model.  Each  of 
these  situations  may  be  coupled  with  a  high-resolution  rendition  from  the  engagement- 
level  model.  Next,  each  situation  in  the  campaign  model  is  linked  to  one  of  the  input 
cases.  The  mean  outcome  of  THUNDER  engagements  is  calibrated  to  approximate  the 
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mean  outcome  of  similar  engagements  in  BRAWLER.  A  scaling  process  allows  any 
engagement  whose  outcome  has  not  been  decided,  or  adjudicated,  in  the  high-resolution 
model  to  be  derived  from  another  engagement  that  was.  THUNDER  uses  a  sequential 
scaling  process  that  takes  combat  between  individual  aircraft  types  and  aggregates  those 
results  up  to  the  flight  group  versus  flight  group  level  for  an  overall  engagement  loss 
score  (AFSAA,  1996).  This  process  is  discussed  in  Section  3.4.4.  Finally,  the  campaign 
model  uses  data  from  the  engagement-level  model  in  order  to  estimate  attrition. 
(Hillestad  and  Moore,  1996) 

3.4.  Combat  Models 

3.4.1.  Introduction 

In  this  section,  we  provide  information  about  the  combat  models  used  in  this 
research.  An  overview  of  BRAWLER  and  THUNDER  is  provided  in  Sections  3.4.1  and 

3.4.2,  respectively.  THUNDER’S  air-to-air  submodel  is  described  in  Section  3.4.4. 
Finally,  an  analysis  of  the  options  considered  by  ASC/XRA  for  calibrating  THUNDER 
with  BRAWLER  is  presented  in  Section  3.4.5. 

3.4.2.  BRAWLER  Overview 

BRAWLER  is  a  high-resolution  air-to-air  combat  simulation  model  that  is  used 
for  engagement-level  analysis  of  few-on-few  air  combat.  It  is  an  event-driven,  Monte 
Carlo  simulation  of  within- visual-range  and  beyond-visual-range  air-to-air  combat 
engagements  between  multiple  flights  of  aircraft,  which  are  explicitly  modeled  through 
all  phases  of  air-to-air  combat.  BRAWLER  is  considered  a  high-resolution  model  due  to 
the  engineering-level  models  of  hardware  and  physical  effects  that  it  incorporates. 
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Typical  uses  of  BRAWLER  include  aircraft  performance  trade  studies,  avionics  and 
weapons  effectiveness  studies,  and  tactics  development  (ASC/XRA,  1996).  BRAWLER 
is  predominantly  written  in  FORTRAN,  with  specific  routines  written  in  C. 

Each  run  of  BRAWLER  simulates  a  single  air-to-air  combat  engagement.  Due  to 
the  stochastic  nature  of  BRAWLER  events,  multiple  replications  of  BRAWLER  using  the 
same  initial  starting  conditions  may  yield  drastically  different  results.  BRAWLER 
produces  results  similar  to  more  detailed  exercises  like  flight  tests  and  man-in-the-loop 
simulations.  Figure  3.2  gives  an  overview  of  BRAWLER. 


Figure  3.2.  BRAWLER  Overview  (ASC/XRA,  1996) 
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3.4.3.  THUNDER  Overview 


THUNDER  is  a  two-sided  theater-level  stochastic  simulation  model  of 
conventional  air,  land,  and  naval  air  warfare.  THUNDER  is  written  in  the  SIMSCRIPT 
II.  5™  programming  language  and  is  used  by  a  large  number  of  US  and  allied  defense 
organizations  and  contractors.  It  is  primarily  used  in  order  to  evaluate  force  structures, 
conduct  Analyses  of  Alternatives  (AOA's),  and  develop  strategies  and  tactics. 

THUNDER  models  air  warfare  at  a  high  level  and  it  also  models  major  aspects  of 
conventional  land  warfare.  Because  it  models  so  many  things  at  a  high  level  of 
aggregation,  THUNDER  is  considered  a  low-resolution  model.  Since  THUNDER  is 
designed  so  that  the  data  is  not  built  into  the  model,  a  THUNDER  baseline  input  database 
is  created  for  a  specific  theater  combat  scenario.  Model  data  is  then  stored  in  a  database 
made  up  of  many  data  files. 

3.4.4.  THUNDER’S  Air-to- Air  Submodel 

THUNDER  contains  its  own  air-to-air  combat  submodel  that  models  air-to-air 
engagements  at  a  very  aggregated  level.  THUNDER  can  be  run  using  a  high-  or  low- 
resolution  setting  for  air-to-air  combat,  which  determines  which  part  of  its  air-to-air 
engagement  submodel  it  uses  to  adjudicate  air  combat.  THUNDER’S  low-resolution  air- 
to-air  setting  allows  only  one  aggregate  PK  to  be  set  for  all  attacking  aircraft/weapon 
versus  target  aircraft  combinations  for  both  Blue  and  Red.  This  limitation  of  the  low- 
resolution  setting  is  unrealistic  and  directly  conflicts  with  the  goal  of  the  calibration 
process.  This  is  because  the  calibration  process  attempts  to  assign  the  correct  PK  values 
to  key  THUNDER  air-to-air  interactions,  which  often  involve  many  different  aircraft 
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from  both  the  Blue  and  Red  sides.  If  the  low-resolution  setting  is  used,  there  is  no  need 
to  perform  the  calibration.  Therefore,  THUNDER’S  high-resolution  setting  is  required 
for  the  air-to-air  calibration  process. 

It  should  be  noted,  however,  that  “high-resolution”  in  the  context  of  THUNDER’S 
air-to-air  model  is  a  relative  term.  THUNDER’S  high-resolution  air-to-air  setting  does 
not  allow  THUNDER  air-to-air  engagements  to  be  adjudicated  at  a  level  of  detail  equal  to 
that  in  BRAWLER.  As  stated  in  Section  1.3,  THUNDER  models  many  aspects  of 
warfare,  so  in  order  to  avoid  being  prohibitively  large,  its  air-to-air  model  is  less  detailed 
than  BRAWLER’S. 

THUNDER’S  high-resolution  air-to-air  combat  process  is  controlled  by  four 
primary  inputs  that  are  listed  and  defined  below. 

1 .  Air-to- Air  PK  —  The  probability  that  a  particular  attacking  aircraft  with  a 
particular  weapon  will  kill  a  specific  target  aircraft  given  that  a  weapon  is 
fired  during  an  encounter  between  the  two.  This  can  be  written  as 

PK(Aircraft+  Weapon  v  Target )  =  P(Killing  Target\Firing  at  Target)  (3.1) 

2.  Missile  Launches  Per  Engagement  --  The  number  of  missile  launches  per 
engagement  for  a  specific  aircraft  weapon  load  for  each  engaged  aircraft. 

3.  Relative  Range  Advantage  —  A  weapon  range  advantage  comparison  between 
opposing  aircraft/weapon  combinations  given  an  encounter  between  the  two 
has  taken  place. 

4.  Detection  Capability  —  The  probability  that  an  aircraft  will  complete  an  air-to- 
air  engagement  against  an  opponent  given  the  opponent  has  been  detected 
(ASC/XRA,  1996). 

ASC/XRA  has  determined  that  inputs  1  and  2  above  may  be  calibrated  using 
BRAWLER.  As  previously  mentioned,  we  are  focusing  on  input  1,  air-to-air  PK. 

In  THUNDER,  an  attacking  aircraft  with  a  specific  weapon  is  called  a  Killer  and  a 
targeted  aircraft  (regardless  of  weapon  type)  is  called  a  Target.  It  is  important  to  note  that 
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THUNDER  pairs  an  attacking  aircraft  with  its  weapon  in  this  manner.  A  PK  value  is 
assigned  to  each  Killer  versus  each  Target.  THUNDER’S  high-resolution  air-to-air 
submodel  adjudicates  air  combat  using  three  steps,  as  follows: 

1 .  A  single  shot  probability  of  kill  (SSPK)  is  calculated  for  each  Killer  versus 
Target  combination.  SSPK  is  defined  for  a  Killer  versus  Target  combination 
as  follows: 

SSPK(Killer  v  Target )  =  PL(Killer)'  PK(Weapon  v  Target )  (3.2) 

where  PL(Killer)  is  the  probability  of  a  weapon  launch  by  the  killer  and 
PK(Weapon  v  Target)  is  the  probability  of  kill  for  the  killer’s  weapon  versus 
the  target  aircraft. 

2.  Next,  single  aircraft  versus  aircraft  SSPK  values  are  aggregated  into  attrition 
rates  for  flight  versus  flight  and  flight  group  versus  flight  group  levels. 

3.  Finally,  engagement  losses  are  assessed  by  random  draws  from  binomial 
distributions. 

3.4.5.  Calibration  of  THUNDER  Using  BRAWLER 

The  process  selected  for  calibrating  THUNDER  using  BRAWLER  is  the  result  of  several 
iterations  of  evaluating  THUNDER  and  BRAWLER  results.  ASC/XRA  ran  several 
experiments  using  THUNDER  and  identified  four  primary  inputs  that  determine  the 
outcome  of  THUNDER  air-to-air  engagements.  They  determined  that  it  is  possible  to 
correlate  two  of  these  four  inputs  with  BRAWLER  results,  as  shown  in  Figure  3.3. 

This  resulted  in  development  of  the  methodology  for  the  calibration  of  THUNDER  with 
BRAWLER.  (ASC/XRA,  1996).  The  process  focuses  on  the  outcome  of  THUNDER 
air-to-air  engagements,  not  whether  or  not  an  engagement  occurs  (Taranto,  1995).  Using 
the  current  THUNDER  calibration  process,  a  BRAWLER  AKPF  MOE  is  used  as  a 
THUNDER  air-to-air  PK  input  value  for  a  Killer  versus  Target  pair.  When  the  Red  and 
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Blue  aircraft  are  the  same  in  each  model  and  the  THUNDER  Killer  fires  the  same  missile 
as  its  BRAWLER  counterpart,  then  this  relationship  can  be  written  as  follows: 

BRAWLER  E(Missile  AKPF)  =  THUNDER  PK  {Killer  v  Target)  (3.3) 

The  AKPF  MOE  is  calculated  in  a  BRAWLER  post-processing  routine  in  a 
procedure  explained  in  Section  4.4.1.  The  BRAWLER  AKPF  MOE  is  the  average 
effectiveness  of  a  missile  including  the  effects  of  missile  kinematics,  tactics,  and  target 
maneuvering  (ASC/XRA,  1996). 


BRAWLER  OUTPUT  THUNDER  INPUT 

<=> 

C=> 


No  BRAWLER  Output 


AIR-TO-AIR  PK 


TOTAL  MISSILE 
LAUNCHES 
PER  ENGAGEMENT 


RELATIVE  RANGE 
ADVANTAGE 


DETECTION  CAPABILITY 


MISSILE 

ADJUSTED  KILLS 
PER  FIRING 


AVERAGE 
MISSILE  LAUNCHES 
PER  ENGAGEMENT 


No  BRAWLER  Output 


Figure  3.3.  BRAWLER  Output  /  THUNDER  Input  Linkages  (Taranto,  1995) 
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3.4.6.  Analysis  of  Calibration  Options 


Two  calibration  options  were  identified  by  ASC/XRA  in  order  to  calibrate  highly 
aggregated  (low-resolution)  THUNDER  air-to-air  engagements  with  high-resolution 
BRAWLER  engagement  analysis,  as  follows: 

1 .  Calibrate  attrition  for  each  distinct  engagement  in  THUNDER,  or 

2.  Calibrate  by  individual  attacker/weapon  performance  versus  an  opponent. 

Option  1  allows  a  direct  correlation  between  BRAWLER  and  THUNDER 

engagements;  however,  the  large  number  of  engagement  types  generated  by  THUNDER 
during  a  multi-day  campaign  makes  running  BRAWLER  for  all  of  these  possible 
engagements  impractical,  perhaps  even  infeasible  (Taranto,  1995).  Additionally,  the 
THUNDER  database  structure  does  not  allow  multiple  PK’s  for  a  Killer  versus  Target 
pair.  This  feature  would  be  required  to  execute  option  1. 

On  the  other  hand,  the  current  THUNDER  database  structure  supports  option  2. 
The  THUNDER  database  is  built  for  aircraft/weapon  versus  aircraft  combinations,  which 
makes  it  necessary  to  calibrate  PK’s  by  individual  missile  system  performance  rather  than 
by  engagement  attrition  (ASC/XRA,  1996).  Option  2  is  far  easier  to  implement,  and  still 
allows  BRAWLER  outputs  and  THUNDER  inputs  to  be  connected,  albeit  in  an  aggregate 
manner. 

Option  2  describes  the  current  THUNDER  calibration  process.  Taranto  (1995) 
documents  several  assumptions  and  limitations  of  the  BRAWLER  calibration  process,  as 
follows: 

1.  A  few  BRAWLER  scenarios  are  assumed  to  be  representative  of  all 
THUNDER  air-to-air  engagements. 
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2.  The  average  effectiveness  of  a  Killer  versus  a  Target  is  approximately  the 
same  regardless  of  whether  the  weapon  is  used  as  the  primary  munition,  used 
with  different  tactics,  or  faces  the  same  opponent  armed  with  different 
weapons. 

3.  The  THUNDER  database  structure  limits  the  air-to-air  PK  inputs  for  an  Killer 
versus  a  Target  to  a  single  value.  For  example,  for  a  MIG-29  armed  with  an 
AA-10A  versus  an  F-15E,  no  corrections  can  be  made  to  the  AA-10A 
effectiveness  when  considering  the  F-15E  armament. 

Based  upon  findings  in  discussions  with  the  THUNDER  and  BRAWLER  analysts 
in  ASC/XRA,  we  should  add  another  limitation  to  this  list,  as  follows: 

4.  BRAWLER’S  calibration  of  THUNDER  is  only  valid  for  the  specific  scenario 
of  interest,  due  to  the  fact  that  the  distribution  of  engagements  may  change 
from  one  scenario  to  another. 

Taranto  alludes  to  this  when  he  states  that  “types  of  engagements  may  change  with 
changes  in  a  concept  of  operations  (CONOPS).”  A  CONOPS  may  be  directly  related  to  a 
scenario  of  interest. 

THUNDER  has  been  calibrated  using  BRAWLER  by  ASC/XRA  most  recently 
using  BRAWLER  Version  6.15.  Individual  air-to-air  engagements  that  occurred  in 
THUNDER  were  compared  to  BRAWLER  MOE’s  and  MOO’s  and  examined  from  two 
different  perspectives.  From  a  THUNDER  input  data  perspective,  the  THUNDER  PK 
input  was  compared  to  BRAWLER  AKPF  output  to  ensure  they  were  similar;  if  not,  an 
appropriate  investigation  was  performed  and  subsequent  corrections  were  made.  From  an 
outcome  perspective,  they  compared  the  outcomes  from  specific  THUNDER 
engagements  to  BRAWLER  engagements  with  similar  force  structures  to  ensure  that  they 
were  also  similar.  In  this  manner,  the  air-to-air  calibration  process  ensures  that 
THUNDER  air-to-air  inputs  and  outputs  are  realistic.  A  simplified  example  of  this 
process  that  focuses  only  on  the  THUNDER  PK  input  is  illustrated  in  Figure  3.4. 
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Figure  3.4.  Calibration  of  THUNDER  Outcomes  and  Input 

A  representative  air-to-air  engagement  force  structure  must  be  chosen  for  each 
Killer/Target  combination.  In  other  words,  one  force  structure,  for  example  4v4,  must  be 
chosen  as  a  baseline  that  all  THUNDER  air-to-air  engagements  between  the  Killer  and 
the  Target  are  derived  from,  using  a  scaling  process.  This  is  necessary  since 
THUNDER’S  database  structure  only  allows  one  PK  for  each  Killer/Target  combination, 
as  previously  mentioned  in  this  section.  This  is  accomplished  by  using  the  most  common 
force  structure  of  an  engagement  for  a  Killer/Target  combination  that  occurs  in 
THUNDER  to  represent  all  expected  engagements  involving  that  Killer/Target 
combination.  The  BRAWLER  and  THUNDER  users  in  ASC/XRA  have  used  this  method 
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in  the  past.  The  BRAWLER  users  polled  the  THUNDER  users  for  the  most  common 
engagement  observed  in  the  THUNDER  scenario  of  interest  for  each  Killer/Target 
combination  (Williams,  1997).  It  should  be  noted  that  the  most  common  force  structure 
of  an  engagement  is  dependent  upon  the  minimum  flight  sizes  that  THUNDER  allows  for 
both  the  Killer  and  Target  aircraft.  A  THUNDER  user  sets  minimum  flight  sizes  by 
aircraft  in  typeac.dat  and  by  mission  in  airplan.dat. 

Recall  from  Section  3.4.4  that  single  aircraft  versus  aircraft  SSPK  values  are 
calculated  for  each  killer  versus  target  combination  and  then  aggregated  into  attrition 
rates  for  flight  versus  flight  and  flight  group  versus  flight  group  levels.  Using  calibration 
option  2,  the  THUNDER  single  aircraft  versus  aircraft  SSPK  values  are  derived  from  the 
individual  Killer  performance  versus  a  Target  in  the  appropriate  representative 
THUNDER  engagement. 

In  the  example  depicted  in  Figure  3.4,  a  4v4  engagement  is  THUNDER’S 
representative  engagement.  Therefore,  the  THUNDER  PK  calibration  using  BRAWLER 
AKPF  is  performed  on  THUNDER  4v4  engagements  using  BRAWLER  4v4 
engagements.  Also,  4v4  BRAWLER  outcomes  are  compared  to  4v4  THUNDER 
outcomes  and  are  directly  calibrated  from  this  comparison. 

4v4  BRAWLER  engagements  may  also  be  compared  to  2v2  THUNDER 
engagements  since  the  force  ratio  is  1:1  for  both  of  these  engagements.  Often,  this  type 
of  comparison  is  made  instead  of  actually  running  the  additional  engagement  in 
BRAWLER  and  making  a  direct  comparison  to  the  similar  engagement  in  THUNDER. 
This  is  because  similar  outcomes  should  be  expected,  and  have  indeed  been  observed  by 
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ASC/XRA,  for  both  of  these  engagements  (Logan,  1998).  A  4v4  BRAWLER 
engagement  and  a  2v2  THUNDER  engagement  are  indirectly  calibrated  in  this  manner. 

In  the  calibration  process,  engagements  in  addition  to  the  representative  THUNDER 
engagement  are  performed  in  BRAWLER,  such  as  2v4.  The  outcomes  of  these 
engagements  may  be  directly  calibrated  as  well.  The  engagements  that  are  directly 
calibrated  are  generally  the  most  common  engagements  in  THUNDER. 
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4.  Methodology 


4.1.  Introduction 

Chapter  4  discusses  the  methodology  of  the  research  in  this  thesis.  Section  4.2  and 
4.3  discuss  specific  details  about  BRAWLER  and  THUNDER,  respectively,  that  apply  to 
this  methodology.  Section  4.4  describes  the  ASC/XRA  THUNDER  calibration  process 
that  was  followed  in  this  research  to  ensure  that  the  THUNDER  PK  input  data  resulting 
from  this  research  has  a  pedigree.  Section  4.5  describes  the  bootstrap  sampling  method 
that  was  used  to  construct  an  alternative  confidence  interval  to  the  standard  confidence 
interval  that  estimates  the  expected  value  of  AKPF.  Finally,  Section  4.6  describes  the 
sequential  procedures  for  constructing  these  two  intervals,  which  we  have  denoted 
Candidate  S 1  and  Candidate  S2. 

4.2.  BRAWLER 

4.2.1.  Introduction 

This  section  highlights  specific  BRAWLER  items  that  were  important  to  this 
research,  beginning  in  Section  4.2.2  with  the  database  employed.  We  also  describe  the 
output  variable  of  interest  (AKPF)  in  Section  4.2.3.  Finally,  we  describe  the  input 
variable  that  was  used  in  this  research  in  Section  4.2.4. 

4.2.2.  Version  and  Database 

Most  BRAWLER  runs  are  made  in  a  classified  environment  when  performing 
analyses  of  weapon  systems.  This  is  due  to  the  classified  nature  of  the  weapon  systems 
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being  modeled.  However,  for  this  research,  an  unclassified  version  of  BRAWLER 
denoted  Version  6.3-  Unclassified  and  an  unclassified  database  provided  by  the  USAF’s 
Survivability/Vulnerability  Information  Analysis  Center  (SURVIAC)  were  used.  The 
unclassified  BRAWLER  database  includes  notional  aircraft,  missiles,  and  radars  for  both 
the  Blue  and  Red  sides.  BRAWLER’S  notional  Blue  aircraft  is  loaded  with  notional 
semi-active  missiles,  similar  to  AIM-7s,  and  notional  infrared  (IR)  missiles  similar  to 
AIM-9s  (Taranto,  1997). 

4.2.3.  Output  Variable 

The  BRAWLER  output  variable  of  interest  is  the  adjusted  kills  per  firing  (AKPF) 
MOE  for  BRAWLER’S  notional  IR  missile.  Kills  per  firing  (KPF)  for  a  missile  type  is 
defined  as  the  proportion  of  the  total  number  of  times  the  missile  killed  its  target  out  of 
the  total  number  of  times  the  missile  was  fired.  To  obtain  AKPF,  we  apply  a  correction 
factor  to  KPF.  Formal  definitions  of  KPF  and  AKPF  are  introduced  in  Section  4.4.1. 

BRAWLER  computes  an  AKPF  for  both  the  IR  and  semi-active  missiles.  We 
decided  to  focus  on  the  AIM-9-like  IR  missile,  since  the  THUNDER  database  used  in  this 
research  contains  the  AIM-9  missile,  but  not  the  AIM-7.  For  the  purpose  of  this  research, 
we  consider  BRAWLER’S  IR  missile  to  be  an  AIM-9.  The  AKPF  value  for  the  AIM-9  is 
used  as  the  THUNDER  PK  input  variable  for  the  F-15C/ AIM-9  aircraft/weapon 
combination,  which  is  discussed  in  Section  4.4.2. 

4.2.4.  Input  Variable 

BRAWLER  contains  three  discrete  pilot  skill  level  settings  that  may  be  used  as 
input  variables:  Rookie  (R),  Pilot  (P),  and  Ace  (A).  Buschor  (1998)  designed  an 
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unclassified  experiment  using  24  different  combinations  of  these  three  pilot  skill  level 
settings  for  the  pilots  of  each  of  four  Blue  aircraft  in  a  4v4  combat  engagement.  In  this 
type  of  engagement,  the  four  aircraft  are  called  a flight,  and  the  aircraft  are  ordered.  This 
particular  size  flight  is  the  ordered  set  of  four  aircraft,  denoted  (1,2, 3,4).  The  skill  level 
combination  of  pilots  in  a  four-aircraft  flight  is  denoted  by  a  set  of  four  letters,  one  for 
each  aircraft  in  the  flight.  For  example,  the  case  where  the  pilots  of  aircraft  1  and  2  are 
Pilots,  aircraft  3  is  an  Ace,  and  aircraft  4  is  a  Rookie  is  denoted  PPAR.  We  call  each 
different  combination  of  input  variables  a  scenario. 

With  Buschor’s  permission,  we  used  BRAWLER  output  data  resulting  from  200 
replications  of  each  of  three  of  his  scenarios,  identified  in  Table  4.1: 

Table  4. 1 .  BRAWLER  Scenarios  Used 


Scenario 

Number 

Pilot  Types 
(1, 2,3,4) 

1 

AAAA 

2 

PPPP 

3 

RRRR 

For  other  BRAWLER  settings  used  in  these  scenarios,  see  Buschor  (1998).  In  this 
research  we  used  this  data  for  the  purpose  of  determining  when  enough  BRAWLER  runs 
have  been  performed  to  obtain  an  accurate  estimate  of  the  mean  AKPF  for  each  scenario. 


4-3 


4.3.  THUNDER 


4.3.1.  Introduction 

In  this  section  we  highlight  the  specific  THUNDER  items  that  were  important  in  this 
research.  First  we  describe  the  THUNDER  database  that  was  used  in  Section  4.3.2.  We 
then  discuss  the  input  variable  of  interest,  the  THUNDER  PK,  in  Section  4.3.3. 

4.3.2.  Version  and  Database 

THUNDER,  like  BRAWLER,  is  run  almost  exclusively  in  a  classified  environment, 
using  classified  input  databases.  A  THUNDER  input  database  is  modified  as  necessary 
to  create  different  scenarios.  During  this  research,  we  used  a  notional,  unclassified 
version  of  a  database  called  STORM  provided  by  ASC/XRA.  STORM  simulates  a 
Southwest  Asia  scenario  similar  to  Operation  Desert  Storm.  The  version  of  THUNDER 
we  used  is  version  6.4.2. 

Because  STORM  is  a  notional  database,  it  is  much  smaller  and  more  generalized 
than  an  actual  THUNDER  database.  For  example,  the  types  of  Blue  Killers  are  limited  to 
“the  major  fighter  aircraft  in  the  US  inventory.  Furthermore,  Blue  Killers  are  grouped 
together  into  a  Blue  Killer  Table  into  classes  depending  on  what  weapon(s)  they  may  fire. 
The  same  is  true  for  the  Red  Killers  and  the  Red  Killer  Table.  Since  the  AIM-9  missile  is 
the  weapon  that  we  are  focusing  on  when  running  BRAWLER,  it  is  appropriate  to  list  the 
Blue  Killers  that  may  fire  an  AIM-9,  which  are  as  follows:  the  F-4G,  EA-6B,  A-10,  F-14, 
F-15C,  F-15E,  F-16,  and  the  F/A-18.  All  of  these  Blue  Killers  with  AIM-9s  are  given  the 
same  killer  identification  number  of  1 1  in  STORM,  which  is  called  a  Killer  ID. 

Similarly,  the  Red  Targets  MIG-29,  MIG-21,  MIG-23,  and  F-l  are  all  assigned  the  same 
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identification  number  of  210,  called  a  Target  ID.  Similar  groupings  of  Red  Killers  and 
Blue  Targets  are  used. 

4.3.3.  Input  Variable 

Air-to-air  PK,  as  defined  in  (3.1),  is  our  THUNDER  input  variable  of  interest.  A  file 
named  airairpkdat  appears  in  the  STORM  database  and  contains  two  air-to-air  PK 
tables.  One  is  designated  the  Blue  PK  Table  and  lists  the  PK  value  for  each  Blue  Killer 
ID  versus  each  Red  Target  ID.  A  similar  table  exists  for  Red  Killers  versus  Blue  Targets, 
designated  the  Red  PK  Table.  In  STORM’s  Blue  PK  Table,  Blue  Killer  1 1  is  assigned  a 
PK  value  of  5  (which  is  shorthand  notation  for  a  PK  of  .50)  versus  Red  Killer  210.  Since 
similar  groupings  of  Blue  killers/Red  Targets  and  Red  killers/Blue  Targets  are  used,  the 
Red  PK  Table  looks  similar  to  the  Blue  PK  table.  STORM’s  airairpkdat  file  appears  in 
Figure  4.1. 

Each  PK  in  the  Blue  and  Red  PK  tables  represents  the  probability  of  a  Killer 
killing  the  Target,  given  a  weapon  firing,  for  all  THUNDER  air-to-air  engagements 
between  the  Killer  and  the  Target.  THUNDER’S  database  structure  only  allows  one  PK 
for  each  Killer/Target  combination,  as  discussed  in  Section  3.4.6.  Also,  the  PK  input 
value  for  the  Killer  versus  the  Target  is  the  same  regardless  of  the  weapon  load  the  Target 
is  carrying  or  the  weapons  the  Target  has  fired. 
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AIR . AIR . PKS .202 

BLUE. PK. MULTIPLIER (DEC)  1.00  RED. PK .MULTIPLIER (DEC)  1.00 

KILLERS 


BLUE 


LOW . RES . AIR . TO . AIR . PKS 
BLUE 

210 

1.000  25 

END. SET 

220 

1.000  25 

END . SET 

230 

1.000  25 

END . SET 

RED 

100 

1.000  2 

END . SET 

110 

1.000  5 

END . SET 

120 

1.000  5 

END . SET 

130 

1.000  5 

END. SET 

END. LOW. RES. PKS 
END. AIR. AIR. PKS 


Figure  4.1.  airairpkdat  (AFSAA,  1995) 
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4.4.  Calibration  of  THUNDER  PK  Using  BRAWLER  AKPF 


4.4. 1 .  Calculating  AKPF 

As  mentioned  in  Section  3.4.6,  the  BRAWLER  output  MOE  of  kills  per  firing 
(KPF)  for  an  engagement  is  transformed  into  the  THUNDER  air-to-air  PK  input  value. 
The  BRAWLER  KPF  output  is  the  average  effectiveness  of  a  missile,  including  such 
effects  as  missile  kinematics,  tactics,  and  target  maneuvering  (ASC/XRA,  1996).  The 
BRAWLER  KPF  output  is  transformed  into  AKPF,  which  becomes  the  THUNDER  PK 
input  for  each  Killer  versus  Target  pair.  This  process  is  accomplished  in  four  steps: 

1 .  A  representative  air-to-air  engagement  force  structure  (i.e.  4v4)  is  chosen  for 
each  Killer/Target  combination. 

2.  A  KPF  MOE  is  obtained. 

3.  An  AKPF  MOE  is  obtained  by  applying  a  correction  factor  to  KPF. 

4.  A  THUNDER  PK  is  assigned. 

In  step  1,  the  force  structure  is  chosen  as  described  in  Section  3.4.6.  Step  2 
requires  running  BRAWLER  a  sufficient  number  of  times  and  calculating  the  mean  KPF 
value  for  those  replications.  The  KPF  MOE  is  calculated  in  a  BRAWLER  post¬ 
processing  routine  created  by  Taranto  called  STATS,  which  appears  in  Appendix  C. 
STATS  processes  each  of  the  BRAWLER  log  files  sequentially  and  totals  all  the  missile 
firings,  fuzings,  and  kills  over  all  the  BRAWLER  replications  performed.  KPF  is  then 
written  to  a  missile  effectiveness  report  in  a  file  called  STATS.prt.  An  example  of  this 
report  appears  in  Appendix  D.  The  KPF  value  is  displayed  as  “Kills/Firing”  in  a  report  in 
STATS.prt  titled  “Missile  Effectiveness  Report.”  The  calculation  of  KPF  is 
straightforward,  according  to: 
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N 

^  Kills 

KPF  =  ^ 1 - 

X  Firings 

'='  (4.1) 

where  N  is  the  number  of  BRAWLER  replications  performed.  Since  all  of  the  equations 

in  this  chapter  involve  sums  calculated  over  N  replications,  we  drop  the  summation 
indices  and  simply  write: 


KPF  = 


^  Kills 
X Firings  . 


(4.2) 


As  shown  in  (4.3),  it  is  not  necessary  to  first  calculate  KPF  to  obtain  AKPF,  but 
calculating  KPF  makes  the  calculation  of  AKPF  more  intuitive. 

In  step  3,  an  AKPF  MOE  is  obtained  by  applying  a  correction  factor  to  KPF.  This 
correction  factor  in  step  3  is  necessary  because  sometimes  a  missile  fuzes  on  a  dead 
target,  alive  when  the  missile  was  fired  but  recently  killed  by  another  missile.  These 
events  are  denoted  FuzingsDEAD  and  are  counted  in  the  STATS  routine.  FuzingsDEAD  are 
not  counted  as  kills  when  calculating  KPF;  however  they  are  counted  as  kills  when 
calculating  AKPF.  The  reason  for  this  is  that  AKPF  measures  the  proportion  of  missiles 
that  fuzed  on  their  targets,  live  or  dead.  Because  of  this,  AKPF  gives  a  better 
representation  of  the  missile’s  overall  effectiveness  than  KPF  does,  and  that  is  what  we 
wish  to  represent  when  we  pass  PK  data  to  THUNDER.  The  formula  for  AKPF  is: 


AKPF  - 


X  Kills 


V  Firing, - 

X Fuzing  s 


(4.3) 
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or  more  intuitively: 


AKPF  =  1 CPF  ■  = - - 

z,  Fuzings  -  2,  FuzingsDEAD  . 

The  derivation  of  (4.4)  from  (4.3)  appears  in  Appendix  H.  Note  that  AKPF  >  KPF 
AKPF  value  appears  in  another  STATS.prt  report  called  “THUNDER  Missile 
Effectiveness  Report.”  As  mentioned  in  Section  3.4.6,  the  AKPF  output  is  only 
representative  of  the  specific  engagement  that  was  run  (i.e.  4v4).  We  are  now  ready  to 
complete  this  process  with  step  4. 

In  step  4,  the  THUNDER  PK  is  assigned  the  resulting  BRAWLER  AKPF  in  step 
3  for  the  THUNDER  representative  engagement.  For  example,  if  the  representative 
engagement  is  4v4,  an  expression  for  the  THUNDER  PK  is  as  follows: 


(4.4) 

The 


PK  =  AKPFm- 


(4.5) 


4.4.2.  Strategy  for  Experiment 

Before  starting  the  air-to-air  calibration  experiment,  we  consulted  a  member  of 
ASC/XRA’s  BRAWLER  analysis  group,  Mr.  Larry  Taranto.  He  suggested  that  in  order 
to  make  our  research  as  realistic  as  possible  with  respect  to  how  the  calibration  process  is 
actually  performed  that  we  take  the  following  steps: 

1.  Investigate  how  many  air-to-air  encounters  the  F-15C  force  typically  makes  in 
the  STORM  scenario  in  THUNDER  to  ensure  that  there  are  a  highly 
significant  number  of  encounters.  If  not,  we  should  consider  choosing  the 
Blue  aircraft  that  had  the  most  air-to-air  encounters.  This  aircraft  is  the  best 
choice  of  a  Blue  aircraft  to  calibrate  first.  This  is  because  the  more  encounters 
an  aircraft  has  determines,  in  part,  its  significance  in  the  campaign.  The  more 
significant  an  aircraft  is  in  the  campaign,  the  more  necessary  it  is  to  calibrate 
this  aircraft’s  THUNDER  engagements  with  BRAWLER  engagements. 
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2.  Determine  which  Red  aircraft  the  F-15C  (if  chosen  first)  most  frequently 
engages.  This  aircraft  is  the  best  choice  of  a  Red  aircraft  to  calibrate  first,  for 
the  same  reason  as  in  step  1  above. 

3.  Let  BRAWLER’S  notional  Blue  aircraft  represent  the  F-15C  and  its  Red 
aircraft  represent  the  F-15C 's  most  frequently  engaged  Red  aircraft.  Then 
data  from  BRAWLER  using  the  notional  weapon  systems  can  be  used  to 
populate  the  notional  F-15C  database  in  THUNDER  (Taranto,  1997). 

Taking  Mr.  Taranto’s  advice,  we  performed  30  THUNDER  replications  of  a  30-day 
war  using  the  STORM  database.  We  determined  that  the  F-15C  has  the  most  air-to-air 
encounters  of  any  Blue  aircraft,  and  it  engages  the  MIG-29  Red  aircraft  most  frequently. 
Therefore,  for  the  purpose  of  this  research,  we  consider  BRAWLER’S  Blue  aircraft  an 
F-15C  and  its  Red  aircraft  a  MIG-29.  Then  BRAWLER  AKPF  data  is  used  to  populate 
the  F-15C/AIM-9  aircraft/weapon  combination  (Killer)  entry  in  the  THUNDER  PK  table. 

4.5.  The  Bootstrap  Sampling  Method 

4.5.1.  Introduction 

Efron  (1979)  popularized  the  idea  of  the  bootstrap  as  a  way  to  estimate  a  sampling 
distribution.  The  bootstrap  is  a  technique  that  uses  computer  processing  power  to 
estimate  the  accuracy  of  statistical  estimates  by  simplifying  certain  complex  calculations 
of  traditional  statistical  theory  (Efron  and  Tibshirani,  1993).  The  bootstrap  method 
involves  analyzing  pseudo-data  obtained  by  repeatedly  and  randomly  resampling  the 
original  data.  In  Section  4.5.2  we  describe  the  bootstrap  sampling  process  and  in  Section 
4.5.3  we  describe  methods  for  constructing  bootstrap  confidence  intervals. 
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4.5.2.  Bootstrap  Sampling 


Much  of  statistical  inference  is  describing  the  relationship  between  a  sample  and 
the  population  from  which  the  sample  was  drawn  (Hall,  1992).  The  bootstrap  method, 
like  other  statistical  inference  methods,  requires  a  random  sample.  A  sample 
X  =  (*j ,  x2 , . . . ,  xn )  is  an  unordered  collection  of  n  data  points  randomly  drawn  from  a 

population  that  can  be  represented  by  a  vector  X  that  represents  the  sample.  The  x,  ’s 

(the  elements  of  the  vector  X )  are  independent  and  identically  distributed  (iid)  random 
variables  with  some  population  distribution  function.  This  population  distribution 
function  of  X  is  denoted  F.  A  sample  that  is  obtained  by  randomly  resampling  with 
replacement  from  the  original  sample  is  called  a  bootstrap  sample.  Formally  defined,  the 
jth  bootstrap  sample  X*  =  (x*,**,. ..,**)  is  an  unordered  collection  of  n  members  of  X 

that  are  drawn  randomly  from  X  with  replacement,  where  the  star  (*)  notation  indicates 
pseudo-data  rather  than  actual  data.  The  fact  that  we  are  sampling  with  replacement 
means  that  X*  may  (and,  in  fact,  almost  certainly  will)  contain  one  or  more  duplicates  of 

the  original  data  points  in  X  .  The  bootstrap  method  can  be  illustrated  with  the  following 
simple  example. 

Suppose  our  original  sample*  =  (0,10,0,40,100,33,33,0,100,50)  represents 
AKPF  - 100%  outputs  from  ten  BRAWLER  runs.  If  we  take  ten  bootstrap  samples 
from  *  ,  we  might  obtain  the  following  results  shown  in  Table  4.2: 


4-11 


Table  4.2.  Ten  Bootstrap  Samples 


. 

* 

* 

* 

* 

* 

* 

* 

* 

* 

J 

_*1 _ 

*2 

*3 

*4 

*5 

*6 

*7 

*8 

Xg 

*10 

1 

0 

33 

100 

100 

40 

100 

40 

33 

40 

50 

2 

40 

10 

0 

33 

10 

100 

50 

0 

0 

100 

3 

33 

33 

100 

0 

33 

50 

40 

0 

0 

10 

4 

40 

10 

40 

33 

0 

33 

0 

10 

0 

0 

5 

33 

10 

50 

100 

33 

50 

100 

50 

50 

33 

6 

40 

33 

33 

50 

10 

100 

0 

50 

100 

50 

7 

33 

50 

0 

0 

100 

40 

0 

100 

33 

100 

8 

100 

33 

50 

100 

0 

0 

50 

33 

33 

33 

9 

33 

0 

100 

0 

0 

10 

100 

33 

33 

10 

10 

33 

33 

100 

50 

0 

33 

40 

100 

100 

0 

When  we  form  bootstrap  samples  in  this  manner,  it  should  be  noted  that  we  are 
treating  the  original  sample  X  as  if  it  represents  a  population.  In  fact,  one  interpretation 
of  the  main  idea  behind  the  bootstrap  is  to  “let  the  sample  play  the  role  of  the  population  ” 
(Kim,  1996).  Another  interpretation  is  that  when  we  use  the  bootstrap  method,  “it  is  as  if 
we  are  sampling  from  an  infinite  population  with  a  composition  that  exactly  matches  that 
of  the  sample  that  we  drew”  (Friedman  &  Friedman,  1995).  So,  when  we  use  the 
bootstrap  method,  we  treat  X  as  if  it  were  the  population  from  which  X  *  was  drawn. 
While  it  may  be  possible  to  compute  the  resampling  distribution  of  X  analytically,  it  is 
often  quicker  and  easier  to  simply  “build  up  enough  of  a  picture”  of  the  distribution  of  X 
by  using  a  simple  computer  simulation  (Ripley,  1987). 

In  a  similar  manner  as  the  population  distribution  function  was  denoted  F,  the 
distribution  function  of  the  sample  X  *  is  denoted  F  .  F  is  called  the  empirical 
distribution  function,  which  is  conditional  on  F.  We  use  the  hat  (A)  notation  to  remind  us 
that  F  gives  us  an  estimate  of  F. 
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The  bootstrap  method  is  often  used  to  estimate  a  parameter,  such  as  the  mean,  of  a 
population  with  a  specific  distribution.  A  parameter  9  of  a  population  is  a  functional  of 
the  distribution  function  F,  which  can  be  denoted  6  =  9(F) .  When  we  treat  our  sample 
like  a  population,  an  estimated  parameter  of  the  distribution  function  of  the  sample 
should  be  a  good  estimator  of  the  parameter  of  the  population.  We  obtain  this  estimated 
parameter  as  follows: 

1 .  Form  a  bootstrap  sample  by  randomly  sampling  n  times,  with  replacement, 
from  the  original  sample  of  size  n. 

2.  Repeat  this  process  B  times,  where  £  is  a  large  number,  so  that  we  obtain  B 
bootstrap  samples  of  size  n. 

3.  For  j  =  1  to  B,  calculate  a  parameter  9 j  of  the  jth  bootstrap  sample  that  is  based 

A  A  A 

on  the  bootstrap  distribution  of  9  =  9(F).  These  parameters  9 j  give  us  a 
good  idea  of  what  the  distribution  of  9  looks  like. 

It  should  be  noted  that  the  bootstrap  method  is  a  stochastic  process  and  if  the 
procedure  is  repeated,  slightly  different  results  should  be  expected.  Let  us  return  to  our 
example  and  suppose  that  we  are  interested  in  estimating  9  ,  where  9  is  the  population 
mean.  Then  we  would  be  interested  in  finding  9j ,  the  mean  of  bootstrap  sample  j  . 

These  values  have  been  added  to  our  original  table  (Table  4.2)  and  appear  in  Table  4.3. 
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Table  4.3.  Ten  Bootstrap  Samples  and  Their  Respective  Means 


n 

* 

*i 

* 

*2 

* 

*3 

*4 

* 

x5 

*6 

* 

*7 

*8 

* 

*10 

A 

0 

33 

100 

100 

40 

100 

40 

33 

40 

50 

53.6 

40 

10 

0 

33 

10 

100 

50 

0 

0 

100 

34.3 

33 

33 

100 

0 

33 

50 

40 

0 

0 

10 

29.9 

40 

10 

40 

33 

0 

33 

0 

10 

0 

0 

16.6 

5 

33 

10 

50 

100 

33 

50 

100 

50 

50 

33 

50.9 

6 

40 

33 

33 

50 

10 

100 

0 

50 

100 

50 

46.6 

33 

50 

0 

0 

100 

40 

0 

100 

33 

100 

45.6 

100 

33 

50 

100 

0 

0 

50 

33 

33 

33 

43.2 

33 

0 

100 

0 

0 

10 

100 

33 

33 

10 

31.9 

33 

33 

100 

50 

0 

33 

40 

100 

100 

0 

48.9 

We  put  our  example  aside  for  now,  but  will  return  to  it  in  Section  4.5.3.4 

The  reader  may  have  reservations  about  this  process,  as  we  did  when  first  faced 
with  the  idea  of  the  bootstrap.  In  fact,  Ripley  (1987)  goes  so  far  as  to  write,  “At  first  sight 
this  is  a  preposterous  idea,  but  it  has  been  shown  to  work  well  in  a  wide  variety  of 
problems.”  The  key  assumption  in  the  bootstrap  method  is  that  the  sample  represents  the 
population.  This  heavy  leverage  on  the  original  sample  is  a  shortcoming  of  the  bootstrap 
method  (Friedman  and  Friedman,  1995).  The  advantage  of  using  bootstrap  sampling  is 
that  we  can  obtain  an  estimate  of  the  parameter  of  a  distribution  of  interest  without 
having  to  make  any  assumptions  about  its  distribution. 

4.5.3.  Bootstrap  Confidence  Intervals 

4.5.3. 1.  Introduction 

Many  statistical  inference  procedures  involve  estimating  a  parameter(s)  9  of  a 
distribution.  An  interval  estimator,  more  commonly  known  as  an  approximate  confidence 
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interval,  is  a  procedure  that  specifies  a  method  for  using  a  sample  to  form  a  specific 
interval.  This  interval  is  designed  such  that,  ideally,  it  possesses  two  attributes: 

1 .  The  interval  contains  the  parameter  9  ,  and 

2.  The  interval  is  relatively  narrow. 

The  endpoints  of  the  interval  are  called  the  upper  and  lower  confidence  limits,  and  since 
they  are  functions  of  a  random  sample,  they  themselves  are  random  quantities.  Therefore, 

we  cannot  be  certain  that  one  interval,  which  is  calculated  from  one  sample,  contains  9  . 

\ 

However,  if  we  use  a  proper  interval  estimator,  we  can  be  confident  that,  when  we 
repeatedly  collect  samples,  the  interval  generated  from  each  sample  will  contain  9  with 
probability  C.  A  confidence  coefficient  C  is  defined  as  the  approximate  fraction  of  the 
time  that,  in  repeated  sampling,  these  intervals  contain  9  ;  or  in  other  words,  the 
approximate  probability  that  the  confidence  interval  will  contain  9  .  Thus,  if  C  for  our 
interval  is  high  (close  to  1),  then  we  can  be  highly  confident  that  a  single  confidence 
interval  calculated  from  a  single  sample  contains  9  (Wackerly,  Mendenhall,  and 
Scheaffer,  1996). 

4.5. 3.2.  Standard  Confidence  Intervals 

The  favorite  approximate  confidence  interval,  by  far,  is  the  standard  confidence 
interval  (DiCiccio  and  Efron,  1996).  Suppose  we  have  one  random  sample  X  from  a 

xv  XV 

population  with  an  unknown  distribution  F.  Furthermore,  let  9  =  9(F)  be  the  estimate  of 

A 

a  parameter  9  =  9(F) ,  and  let  (ffl-  be  the  standard  error  of  9  .  The  standard  error  is  a 
useful  measure  of  the  accuracy  of  9  as  an  estimate  of  9  .  Standard  errors  are  commonly 
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used  to  develop  approximate  confidence  intervals  for  9  .  An  expression  for  d  .  is 

6 

developed  in  the  equations  that  follow. 

The  standard  deviation  6  of  our  sample  X  of  xt-  ’s  is  defined  as 


d  = 


(4.6) 


where  x  is  the  mean  of  the  sample  and  n  is  the  sample  size.  We  continue  to  use  the  hat 
notation  to  indicate  that  d  is  an  estimate  of  the  population’s  standard  deviation  c .  From 
this  point  forward,  we  shall  assume  that  the  parameter  9  that  we  are  trying  to  estimate  is 

/V 

the  mean,  so  in  this  case,  9  =  x .  Then  we  can  replace  (4.6)  with  the  following: 


<7  = . 


]n-\ tf 


■9f 


(4.7) 


The  standard  error  of  9  is  equivalent  to  the  standard  deviation  of  9  ,  which  we  have 
previously  denoted  dg- .  An  expression  for  d-  is 


(4.8) 


The  so-called  large-sample  result  maintains  that  under  most  circumstances,  as  the 
sample  size  n  grows  large,  the  distribution  of  9  approaches  normality,  with  mean  near  9 
and  variance  near  d 2 ,  written  9  ~  N(9,  d 2 ) .  From  this  result  we  define  Z  as  the 
following: 


Z  = 


9-9 


(4.9) 
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where  Z  ~  N(0,1)  as  the  sample  size  We  say  that  Z  is  distributed  according  to  the 

standard  normal  distribution  or  z  distribution.  The  standard  confidence  interval  is 
derived  from  the  assumption  that  Z  ~  N(0,1).  However,  this  assumption  is  only  valid  for 


an  infinite  sample  size.  For  6  =  x ,  which  we  have  already  assumed  to  be  the  case,  we 
use  a  better  approximation  derived  by  Gosset  (1908): 


(n-1), 


(4.10) 


which  is  distributed  according  to  the  Student’s  t-distribution  with  n-1  degrees  of  freedom. 
We  may  now  define  a  100  •  (1  -a)%  standard  confidence  interval  as  follows: 


6±t, 


(n-1, a/2) 


or  equivalently  as 


(4.11) 


t(n-l,a/2) 


■6s,6+tt 


(n-1, a/2) 


(4.12) 


where  a  =  1  -C  is  called  the  level  of  significance  and  f(n_la/2)  is  the  100  a/2  percentile 


of  a  t  distribution  with  n-1  degrees  of  freedom.  We  call  the  distance  from  0  to  each 
endpoint  of  the  interval  the  half-length  of  the  confidence  interval  and  denote  it  as 


hl(n,oc)  —  <5$ 


(4.13) 


where  hi  is  a  function  of  n,  the  sample  size,  and  a ,  the  significance  level. 

To  help  transition  towards  how  (and  why)  the  bootstrap  method  is  used  to 
generate  confidence  intervals,  let  us  consider  representing  the  100  •  (1  -  a)%  interval  as 
follows: 
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[0(a/2)  0(l~al2) 


(4.14) 


where  0ia/2->  and  0 11  a/2'  are  the  1 00 •  cnj2  and  100-(1  —  oc/2)  percentiles,  respectively, 

of  the  distribution  of  9  .  Now  we  may  define  the  length  of  a  100  •  (l-a)%  confidence 
interval  as 


length  =  9[l-  a  /  2]-  9[a /2] 


(4.15) 


We  may  also  define  the  shape  of  the  interval  in  terms  of  asymmetry  around  9  ,  or 
alternatively,  skewness,  as 


shape = 


0[l-a/2]-g 
9-9[al  2] 


(4.16) 


Therefore,  an  interval  that  is  more  dense  to  the  left  of  9  (skewed  right)  has  a  shape 
greater  than  1,  while  the  opposite  case  (skewed  left)  has  a  shape  less  than  1.  A  symmetric 
interval  has  a  shape  equal  to  1.  We  have  already  seen  from  (4.1 1)  that  the  standard 
confidence  interval  is  symmetric. 

Standard  confidence  intervals  always  have  shape  equal  to  1,  and  it  is  in  this  way 
that  they  can  be  quite  inaccurate  in  practice  (DiCiccio  and  Efron,  1996).  Using  a  standard 
confidence  interval  requires  that  we  make  assumptions  based  on  a  normal  distribution, 
which  may  not  be  appropriate.  We  show  in  Sections  4.5.3.3  and  4.5.3.4  that  this 
symmetric  property  is  not  necessarily  the  case  with  bootstrap  confidence  intervals.  This 
characteristic  is  a  potential  asset  of  the  bootstrap  method. 
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4.5.3. 3.  Bootstrap-/  Confidence  Intervals 


The  bootstrap  sampling  method  may  be  used  to  construct  an  approximate 

A 

confidence  interval  for  an  estimated  parameter  6 .  This  parameter  is  estimated  by  a 
statistic  that  is  a  function  of  n  independent  and  identically  distributed  (iid)  observations. 
The  process  to  construct  a  bootstrap-?  confidence  interval  begins  as  follows: 

1 .  Generate  B  independent  bootstrap  samples  of  size  n  using  the  empirical 
distribution  of  the  original  n  observations. 

2.  Obtain  the  empirical  distribution  of  the  statistic  from  these  B  independent 
bootstrap  samples. 

For  sufficiently  large  B  one  obtains  an  arbitrarily  good  approximation  to  a  distribution 
that  is  a  good  estimate  of  the  true  distribution  of  the  statistic. 

By  using  the  bootstrap  procedure  we  can  obtain  accurate  intervals  without  having 
to  make  assumptions  about  a  standard  normal  distribution  as  in  (4. 10).  The  idea  of  the 
bootstrap-?  interval  is  to  estimate  the  percentiles  of  Zby  bootstrapping.  The  bootstrap-? 
interval  procedure  estimates  the  distribution  of  Z  by  bootstrap  sampling  from  the  data. 

By  doing  this,  the  bootstrap-?  interval  automatically  adjusts  for  skewness  in  the 
underlying  population.  We  compute  a  bootstrap  value  of  Z  for  each  bootstrap  sample  b 
using  the  following  equation: 


hL 

(&) 


(4.17) 


where  0*(b)  is  the  estimate  of  6  for  the  bootstrap  sample  x*b ,  <t‘-  is  the  estimated 
standard  error  of  6*  for  the  bootstrap  sample  x*b  and  §*=  (l/fi)  •  ^fh_{Q*(b)  (DiCiccio 
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and  Efron,  1996).  The  (a/2-100)  percentile  (%)  of  Z*(b)  is  estimated  by  the  value  ta/2 
such  that 


Z\b)<tan=ct/2-B 

(4.18) 

where  a  is  the  significance  level. 

The  estimate  of  the  (a  /  2  •  100)%  point  of  Z  is  the  (a  /  2)  •  B  ordered  value  of 
Z*(b)  and  the  estimate  of  the  [(1 -a/2) -100)]%  point  is  the  (l-a/2)-F  ordered  value 
of  Z*  (b) .  The  bootstrap-?  confidence  interval  then  is 


K®  *<l-o/2)  )>  (0  *(a/2)'^e)] 

4.5.3.4.  Bootstrap  Percentile  Intervals 


(4.19) 


The  bootstrap  percentile  interval  method  is  another  simple  method  for  assigning 
an  approximate  confidence  interval  to  a  real-valued  parameter  6  =  9(F)  that  is  based 

upon  the  bootstrap  distribution  of  9  =  9(F) .  In  section  4.5.3.2  we  considered 
representing  the  100-  (l-a)%  confidence  interval  in  the  form  of  Equation  (4.15).  Let 


A  A^  ^  ^ 

F(t)  =  P{9  <  ?}  be  the  empirical  distribution  function  of  the  bootstrap  distribution  of  9  . 
Wethen  define  0(“/2)  =F_1(a/2)  and  0(1~“/2)  =F"‘(l-a/2).  The  percentile  method 
allows  us  to  take 


(l-a/2) 


(4.20) 


as  an  approximate  100-(l-a)%  central  confidence  interval  for  9  (Efron,  1982).  This  is 
illustrated  in  Figure  4.2. 
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Median  of 

Bootstrapped  Means 


Figure  4.2.  Bootstrap  Percentile  Interval 

If  B  =  501  and  oc=  .10,  which  is  the  case  in  the  program  enoughruns,  then  0ia,2)  is  the 

26th  ordered  value  of  the  replications  and  0(1~“/2)  is  the  476th.  We  discuss  this  further  in 
Section  4.6.3. 

Let  us  return  to  our  example  and  suppose  that  we  wish  to  find  an  80%  confidence 
interval  for  6  ,  the  population  mean.  We  start  by  sorting  our  10  samples  by  the  value  of 

A 

6j ,  the  mean  of  bootstrap  sample  j,  as  shown  in  Table  4.4. 
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Table  4.4.  Ten  Bootstrap  Samples  Sorted  by  Mean 


fl 

* 

* 

*2 

* 

*3 

* 

*4 

* 

*5 

*6 

* 

*7 

*8 

* 

*9 

*10 

o, 

4 

40 

10 

40 

33 

0 

33 

0 

10 

0 

0 

16.6 

3 

33 

33 

100 

0 

33 

50 

40 

0 

0 

10 

29.9 

9 

33 

0 

100 

0 

0 

10 

100 

33 

33 

10 

31.9 

2 

40 

10 

0 

33 

10 

100 

,50 

0 

0 

100 

34.3 

8 

100 

33 

50 

100 

0 

0 

50 

33 

33 

33 

43.2 

7 

33 

50 

0 

0 

100 

40 

0 

100 

33 

100 

45.6 

6 

40 

33 

33 

50 

10 

100 

0 

50 

100 

50 

46.6 

10 

33 

33 

100 

50 

0 

33 

40 

100 

100 

0 

48.9 

5 

33 

10 

50 

100 

33 

50 

100 

50 

50 

33 

50.9 

1 

0 

33 

100 

100 

40 

100 

40 

33 

40 

50 

53.6 

A 

A  80%  bootstrap  percentile  interval  for  0  gives  us  an  approximate  80% 
confidence  interval  for  6  .  We  can  construct  a  bootstrap  80%  percentile  interval  for  6 
using  the  10th  percentile  and  the  90th  percentile  of  the  0y  ’s,  which  are  16.6  and  50.9, 

respectively.  An  80%  bootstrap  percentile  interval  for  6  based  on  our  10  bootstrap 
samples  is  (16.6,50.9).  However,  in  practice,  we  would  never  use  such  a  small  number  of 
bootstrap  samples,  since  it  is  desirable  (and  very  easy)  to  obtain  a  large  number  B  of 
bootstrap  samples.  In  comparison,  a  standard  80%  confidence  interval  for  6  based  on 
our  original  sample  of  size  10  according  to  (4.1 1)  is  36.6  ±  16.6  =  (20.0,53.2) . 


4.6.  Sequential  Procedures  for  Determining  Run  Length 
4.6.1.  Introduction 

In  this  section  we  discuss  how  sequential  procedures  for  determining  simulation 
run  length  are  applied  to  the  problem  at  hand:  obtaining  a  sufficiently  precise  estimate  of 
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AKPF  to  use  to  calibrate  a  THUNDER  PK.  First,  Section  4.6.2  introduces  a  sequential 
procedure  presented  in  Law  and  Kelton  (1991)  that  makes  use  of  a  standard  confidence 
interval.  Section  4.6.3  introduces  a  sequential  procedure  involving  bootstrapping  that  is 
used  in  Johnson  and  Lucas’  (1994)  program  enoughruns  to  determine  when  enough 
BRAWLER  runs  have  been  performed  to  obtain  an  accurate  estimate  of  exchange  ratio. 
This  procedure  is  a  predecessor  to  the  procedure  used  in  Candidate  S2.  Section  4.6.4 
formally  introduces  Candidate  S2,  which  is  based  on  the  procedure  in  enoughruns. 

Finally,  Section  4.6.5  discusses  how  to  prepare  BRAWLER  data  for  use  with  Candidate 
S2. 

4.6.2.  Candidate  SI 

Law  and  Kelton  (1991)  present  a  sequential  procedure  for  obtaining  an  estimate  of 
the  mean  6  with  a  specified  relative  error  y  and  a  confidence  level  a  that  allows  the 
simulation  user  to  only  perform  as  many  replications  as  necessary.  We  choose  an  initial 
number  of  replications  «0>2  and  perform  the  following  procedure: 

1 .  Perform  n0  replications  of  the  simulation  and  set.  n  =  n0. 

2.  From  the  sample  of  outputs  (x, ,  x2 , . . . ,  xn ) ,  compute  the  sample  mean 
x(n)  and  the  half-length  of  the  confidence  interval  hl(n,a) . 

3.  Determine  if  the  following  inequality  holds: 

hl(n1cO<y, 

|*(")|  ”  ’  (4.21) 

where  e  is  the  observed  relative  error  and  y'  =  y/(y  +  1)  is  the  adjusted  relative  error 

needed  to  obtain  an  actual  relative  error  of  y .  If  so,  stop  performing  replications  and  use 
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x(n)  as  an  estimate  for  6  .  If  not,  set  n  =  n  +  \,  perform  another  replication,  and  return 
to  step  2. 

Law  and  Kelton  recommend  using  this  sequential  procedure  with  n0  >  10  and 
7  <  0.15 .  After  testing  this  procedure  on  many  different  simulation  models,  they  found 
that  if  these  recommendations  are  followed,  the  estimated  proportion  of  time  that  9  was 
covered  using  a  90%  confidence  interval  was  never  less  than  0.864. 

4.6.3.  Candidate  S2  Predecessor 

Johnson  and  Lucas’  program  enoughruns  is  the  predecessor  to  the  procedure  that 
we  will  use  to  construct  a  bootstrap  percentile  interval  for  the  AKPF  MOE.  The  program 
focuses  specifically  on  exchange  ratios  and  automatically  determines  if  enough 
BRAWLER  runs  have  been  completed,  based  on  the  stopping  criterion  defined  later  in 
this  section.  The  program  enoughruns  is  automatically  executed  between  each 
BRAWLER  run  after  20  runs  have  been  performed. 

Recall  that  the  equation  for  calculating  the  exchange  ratio  for  a  replication  i  ( ERi ) 
is  given  by  (2.1).  To  use  bootstrap  sampling,  we  must  be  able  to  form  bootstrap  samples 
from  an  original  sample  of  data  from  individual  runs.  The  ERi  ’s  from  each  replication  i 

are  such  a  sample.  It  might  seem  that  (2.1)  would  allow  us  to  bootstrap  sample  ERt ’s. 
However,  referring  to  (2.1),  it  should  be  clear  that  ERt  has  a  lower  bound  zero  and  no 
upper  bound,  written  ERe  [0,°°) ,  since  as  BKi  — >  0,  ER,  — >  oo . 

In  the  BRAWLER  scenarios  used  in  this  research,  as  well  as  many  other 
BRAWLER  scenarios,  one  or  more  replications  result  in  an  unbounded  ERt .  Since,  in 
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general,  we  can  not  form  a  sample  of  N  bounded  ERi  ’s,  we  should  not  attempt  to 
bootstrap  a  sample  of  ERl  ’s. 

We  observe  that  ERt  is  simply  an  aggregate  value  that  represents  the  two  values 
RKj  and  BKr  We  form  a  sample  of  N  ordered  pairs  (RKi ,  BKt )  where  ie  [1,/V]  and 
denote  this  sample  ( RK,BK ).  We  bootstrap  from  this  sample,  forming  B  bootstrap 
samples.  The  bth  bootstrap  sample  (RK*  (b),  BK*  (£>))contains  the  ordered  pairs 
( RK *  (b),  BK*  (b))  where  b  e  [1,  B]  and  i  e  [1,  N] .  We  then  obtain  the  mean  exchange 
ratio  ER *  (b)  for  each  of  the  b  bootstrap  samples  by  (4.22),  which  is  based  on  (2.2). 

f^RKjib) 

ER\b)= -f - 

(4.22) 

1=1 

This  is  the  technique  that  Johnson  and  Lucas  use  in  enoughruns.  In  their  program, 
the  following  steps  occur: 

1 .  B  =  50 1  bootstrap  samples  (RK'  (b),  BK*  (b))  of  size  N  x  2  are  formed  from 
the  original  sample  of  N  ordered  pairs  (RKi ,  BKj ) . 

2.  The  mean  exchange  ratio  ER*  (b)  for  each  of  the  501  bootstrap  samples  is 
calculated  using  (4.22). 

3.  The  501  ER*  (b)  ’s  are  sorted  in  ascending  order. 

4.  The  following  statistics  are  calculated  for  the  501  ER*  ( b )  ’s. 

c.  The  median,  or  50th  percentile,  ER*  (b)(  50) ,  which  is  the  251st  ordered 
value  of  the  sorted  ER*  (b)  ’s, 

d.  The  5th  percentile,  ER*(b){05) ,  which  is  the  26th  ordered  value,  and 
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e.  The  95th  percentile,  ER*(b)(95) ,  which  is  the  476th  ordered  value. 

5.  Stop  performing  runs  if  both  ER*  (b)(05)  and  ER*(b)(  95)  are  within  ±10%  of 
ER*  ( b ) ( 50) .  Otherwise  continue. 

In  other  words,  this  stopping  criterion  is  met  when  a  90%  bootstrap  percentile  interval 
with  a  relative  error  less  than  .10  is  obtained. 

The  procedure  is  completely  automated.  If  the  stopping  criterion  is  met, 
enoughruns  returns  a  value  of  “yes,”  indicating  that  enough  runs  have  been  performed 
and  BRAWLER  replications  are  terminated.  If  the  stopping  criterion  is  not  met, 
enoughruns  returns  a  value  of  “no,”  and  another  BRAWLER  replication  is  launched. 
Since  the  bootstrap  method  is  a  stochastic  process,  we  may  infer  that  if  this  procedure 
were  repeated  on  the  same  set  of  BRAWLER  data,  a  different  result  (“yes”  versus  “no”) 
may  be  obtained  on  the  same  run.  The  procedure  that  Johnson  and  Lucas  use  in 
enoughruns  is  illustrated  in  Figure  4.3. 
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Figure  4.3.  Procedure  Used  in  enoughruns 
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4.6.4.  Candidate  S2 


The  procedure  applied  in  this  research  is  based  on  Johnson  and  Lucas’  procedure 
in  enoughruns  that  was  discussed  in  Section  4.6.3;  therefore,  the  two  procedures  are  very 
similar.  However,  there  is  a  key  difference  in  the  way  we  form  a  sample  to  bootstrap 
from.  This  is  because  we  must  consider  not  only  the  AKPF  of  each  run,  but  also  the 
number  of  firings  in  each  run  (Firings).  This  can  be  illustrated  with  a  simple  example. 
Suppose  we  have  three  BRAWLER  AKPF  MOE’s  from  three  different  BRAWLER 
outcomes,  and  the  MOO’s  of  interest  are  the  number  of  kills  (Kills)  and  Firings.  This 
data  appears  in  Table  4.5: 

Table  4.5.  MOEs  and  MOOs  from  Three  Different  BRAWLER  Outcomes 


MOE 

MOO 

AKPF 

Kills 

BM1 

0 

0 

l 

0 

0 

10 

0 

0 

0 

Although  all  of  these  outcomes  result  in  an  AKPF  MOE  of  zero,  they  are  clearly  different. 
If  we  only  consider  the  AKPF  values  from  each  run,  we  might  attempt  to  bootstrap 
AKPF’s  from  a  sample  of  size  N  =  3  (the  number  of  BRAWLER  outcomes  or  runs).  If 
we  do  this,  we  give  each  of  these  three  AKPF’s  equal  weight,  which  they  clearly  should 
not  have.  Rather,  an  AKPF  based  on  10  firings  should  carry  10  times  the  weight  of  an 
AKPF  based  on  1  firing,  and  an  AKPF  based  on  zero  firings  should  carry  no  weight. 

We  need  to  collect  two  values  for  each  run,  AKPF  and  Firings,  and  store  them  in 
an  A  x  2  matrix.  We  then  create  one  vector  of  AKPF’s  that  accounts  for  the  number  of 
firings  associated  with  each  AKPF.  This  is  a  key  difference  from  Johnson  and  Lucas’ 
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procedure.  We  use  our  N  x.2  matrix  to  create  an  equal  weight  (with  respect  to  Firings) 
vector  of  AKPF’s  of  size  V*  Firings  .  In  this  vector,  which  is  denoted  the  AKPF 

vector,  the  AKPF  for  each  replication  is  repeated  a  number  of  times  equal  to  the  number 
of  Firings  in  that  replication,  which  allows  us  to  preserve  information  about  the  number 
of  missile  firings  in  each  replication.  We  now  treat  the  new  AKPF  vector  as  our  original 
sample  and  bootstrap  from  this  vector.  For  example,  AKPF  output  from  the  first  five 
replications  of  a  notional  BRAWLER  scenario,  denoted  Scenario  X,  appears  in  Table  4.6. 
For  illustrative  purposes  we  have  also  displayed  the  Kills  in  each  run. 


Table  4.6.  Kills,  Firings  and  AKPF  Matrix 


Replication 

Number 

Scenario  X 

Kills 

F5ff!^3 

AKPF 

1 

1 

2 

0.50 

2 

0.00 

3 

1 

2 

1.00 

4 

1 

4 

0.25 

5 

2 

0.00 

Note  that  in  replication  number  3,  although  there  were  2  firings  but  only  1  kill,  the  AKPF 
is  1.00.  This  is  due  to  the  fact  that  in  this  replication,  one  of  the  missiles  fired  killed  the 
target  and  the  other  fuzed  on  the  same  target  after  it  was  dead. 

We  wish  to  form  one  vector  of  AKPF’s  that  accounts  for  the  Firings  associated 
with  each  AKPF.  As  an  intermediate  step,  Table  4.7  shows  the  results  of  repeating  the 
AKPF  in  each  replication  a  number  of  times  equal  to  Firings  in  that  replication. 
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Table  4.7.  Forming  an  AKPF  Vector 


Replication 

Number 

Scenario  X 

Firings 

AKPF 

1 

1 

0.50 

1 

1 

0.50 

3 

1 

1.00 

3 

1 

1.00 

4 

1 

0.25 

4 

1 

0.25 

4 

1 

0.25 

4 

1 

0.25 

5 

1 

0.00 

5 

1 

0.00 

Note  that  since  replication  2  had  no  Firings,  its  AKPF  of  zero  is  not  entered  into  the  new 
AKPF  vector.  Finally,  we  display  the  vector  of  AKPF’s  that  we  may  bootstrap  from  in 
Table  4.8. 


Table  4.8.  AKPF  Vector  for  Bootstrap  Sampling 


AKPF 

0.50 

0.50 

1.00 

1.00 

0.25 

0.25 

0.25 

0.25 

0.00 

0.00 


The  process  of  constructing  an  AKPF  vector  is  summarized  in  Figure  4.4. 


BRAWLER  Data 


R  e  p  lie  a  tio  n 

N  u  m  b  e  r 

Scenario  X  1 

K  ills 

F  ir in  a  s 

AKPF 

1 

1 

2 

0.50 

2 

0 

0 

0.00 

3 

1 

2 

1  .00 

4 

1 

4 

0  .25 

5 

0 

2 

0  .00 

Form  an  equal  weight  vector 


Replication 

Number 

Seen* 

ario  X 

AKPF 

Firinas 

1 

0.50 

1 

1 

0.50 

1 

3 

1.00 

1 

3 

1.00 

1 

4 

0.25 

1 

4 

0.25 

1 

4 

0.25 

1 

4 

0.25 

1 

5 

0.00 

1 

5 

0.00 

1 

Bootstrap  from  this  vector 


Figure  4.4.  Creating  an  AKPF  Vector  for  Bootstrapping 
As  in  Johnson  and  Lucas’  program,  we  construct  501  bootstrap  samples,  which  are 

n 

all  AKPF  vectors  of  size  ^  Firings ,  ,  which  we  abbreviate  as  ZF.  Once  we  have  501 

i=l 

bootstrap  samples  of  size  ZF,  we  take  the  following  steps,  which  are  similar  to  the  steps 
followed  in  enoughruns. 

1.  Form  501  bootstrap  samples  of  size  ZF  from  the  AKPF  vector. 

2.  Calculate  the  mean  AKPF,  AKPF*(&) ,  of  each  of  the  501  bootstrap  samples. 

3.  Sort  the  501  AKPF*  (b)  ’s  in  ascending  order. 

4.  Calculate  the  following  statistics  for  the  501  AKPF*  (b)  ’s: 

a.  The  median,  or  50th  percentile,  AKPF  *  (b) ( 50) ,  which  is  the  25 1  st  ordered 
value  of  the  sorted  AKPF*  (b)  ’s, 

b.  The  5th  percentile,  AKPF* (b)<05) ,  which  is  the  26th  ordered  value,  and 

c.  The  95th  percentile,  AKPF*  (b)(95) ,  which  is  the  476th  ordered  value. 

5.  Stop  performing  runs  if  both  AKPF*  (b)i05)  and  AKPF*  (b){  95)  are  within 
±10%  of  AKPF’(fc)(50).  Otherwise  continue. 
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4.6.5.  Preparing  BRAWLER  Data  for  Candidate  S2 


In  order  to  apply  Candidate  S2,  consider  how  must  obtain  the  AKPF  from  each 
individual  BRAWLER  run.  We  do  this  by  running  Taranto’s  STATS  post-processing 
program  on  each  individual  BRAWLER  log  file.  However,  STATS  is  designed  to 
calculate  the  mean  AKPF  for  all  BRAWLER  replications  that  have  been  performed,  and 
it  automatically  and  sequentially  processes  all  the  BRAWLER  log  files  that  appear  in  a 
working  directory  and  then  returns  the  mean  AKPF.  So,  in  order  to  obtain  the  AKPF  for 
each  individual  log  file  in  a  directory  instead  of  all  the  runs  considered  en  masse,  we 
created  a  simple  Unix  C-Shell  program  called  make _flyAKPF ,  found  in  Appendix  F. 

This  program  calculates  the  AKPF  values  for  each  run  individually  using  STATS  and 
creates  one  output  file  of  AKPF  and  Firings  from  each  run.  This  file  is  easily  transformed 
into  an  iV  x  2  matrix  for  use  with  S2.  The  process  for  preparing  BRAWLER  AKPF  data 
for  use  by  S2  is  displayed  in  Figure  4.5. 

Conventional  Method: 


Bootstrap  Method: 


Figure  4.5.  Preparing  BRAWLER  AKPF  Data  for  Bootstrapping 
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5.  Results 


5.1.  Introduction 


Recall  that  the  main  purpose  of  this  research  is  to  investigate  if  the  number  of 
BRAWLER  replications  necessary  to  obtain  a  bootstrap  confidence  interval  that  meets 
our  stopping  criterion  is  less  than  the  number  required  to  obtain  a  standard  confidence 
interval  meeting  the  same  criterion.  By  setting  a  =  .10  and  y  = .  10 ,  our  stopping  criteria 
is  as  follows: 

•  Stop  performing  replications  when  a  90%  confidence  interval  about  the 
estimated  mean  can  be  constructed  such  that  both  endpoints  of  the  interval 
are  within  ±  10%  of  the  mean,  or  in  other  words,  the  90%  confidence  interval 
satisfies  the  requirement  that  the  relative  error  of  the  mean  is  .10. 

The  values  a  = .  10  and  y  = .  10  were  chosen  primarily  because  they  are  the  same 
values  that  are  used  in  enoughruns.  In  accordance  with  the  objective  of  the  investigation, 
the  following  steps  were  performed  for  each  of  the  three  scenarios. 

1 .  200  BRAWLER  replications  were  performed  for  each  scenario  and  the  90% 
confidence  intervals  for  the  mean  AKPF  were  calculated.  These  confidence 
intervals  are  our  Truth  model . 

2.  The  minimum  number  A,  of  BRAWLER  runs  necessary  to  provide  a  90% 
standard  confidence  interval  about  the  estimated  mean  with  a  relative  error  of 
.  10  for  the  mean  AKPF  was  found  for  each  scenario  using  Candidate  SI.. 

3.  By  bootstrapping  the  AKPF  data  B  times,  where  B  =  501,  the  minimum 
number  N2  of  BRAWLER  replications  necessary  to  provide  a  bootstrap  90% 

percentile  interval  that  meets  the  same  criterion  as  in  part  2  was  found  for  each 
scenario  using  Candidate  S2. 

4.  N]  and  N2 ,  along  with  the  two  intervals  generated  by  each  method  were 

compared  for  each  replication  of  each  scenario.  The  following  questions  were 
answered: 
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a.  IsE(A2)<  E(  A, )?  In  other  words,  can  a  90%  bootstrap  percentile 
interval  using  N2  replications,  where  N2  <  N]  be  constructed  that  meets 
the  same  stopping  criterion  as  a  standard  90%  confidence  interval  using 
A,  replications? 

b.  Do  either  (or  both)  of  the  two  intervals  calculated  in  steps  2  and  3  overlap 
the  interval  calculated  in  step  1?  In  other  words,  were  the  estimates 
statistically  similar? 

An  exception  to  step  1  was  Scenario  2,  in  which  BRAWLER  aborted  an  excessive 
number  of  replications  for  various  reasons  so  that  only  193  replications  were  performed. 
The  AKPF  and  Firings  data  for  each  run  appears  in  Appendix  G.  The  results  from  step  1 
appear  in  Table  5.1,  where  Cl  Low  and  Cl  Up  are  the  lower  and  upper  bounds  of  the  90% 
confidence  interval,  respectively.  The  results  from  steps  2  and  3  are  presented  in  Sections 
5.3  through  5.5. 


Table  5.1.  Truth  Models  Resulting  from  200  BRAWLER  Runs 


Scenario 

#  Runs 

Cl  Low 

Mean 

Cl  Up 

Cl  Length 

1 

200 

30.07 

32.18 

34.28 

4.21 

2 

193 

29.49 

31.54 

33.59 

4.09 

3 

200 

26.01 

28.49 

30.97 

4.96 

5.2.  Candidate  Measures  of  Merit 

In  this  section  we  discuss  certain  properties  that  we  deemed  appropriate  to  assess 

the  merit  of  each  candidate.  We  call  such  a  property  a  measure  of  merit  (MOM).  The 

two  MOMs  that  we  chose  to  evaluate  Candidate  i,  r  =  1,2  were  as  follows: 

1.  E(lVj):  The  expected  value  of  the  number  of  runs  (AO  required  to  meet  our 
stopping  criterion. 

This  MOM  is  straightforward.  The  smaller  the  value  of  E(A0,  the  better,  since  the  more 
runs  we  must  perform,  the  more  time  and  computer  resources  we  must  use. 
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2.  Probability  of  Overlap  (Pi):  This  probability  is  estimated  by  the  proportion  of 
time  the  confidence  intervals  constructed  using  Candidate  i  and  the  confidence 
interval  constructed  from  our  Truth  model  overlapped. 

When  confidence  intervals  overlap  in  this  way  we  conclude  that  the  results  obtained  from 

using  both  methods  are  statistically  similar.  However,  if  the  intervals  do  not  overlap,  we 

conclude  that  the  results  are  statistically  different. 

Given  these  facts,  we  wish  to  avoid  the  event  that  a  candidate  returns  a  confidence 
interval  that  meets  the  stopping  criterion  but  does  not  overlap  the  confidence  interval 
constructed  from  our  Truth  model .  If  this  event  were  to  occur,  the  interval  returned 
would  be  deceptive  because  if  the  interval  does  not  overlap  the  interval  from  our  Truth 
model ,  it  most  likely  does  not  cover  the  expected  AKPF.  Therefore,  we  define  a  binary 
variable  Overlap.  If  the  interval  constructed  by  each  candidate  for  each  scenario  overlaps 
the  similarly  constructed  interval  for  our  Truth  model ,  we  assign  Overlap  the  value  of 
one,  otherwise,  it  is  assigned  zero  (1  =  overlap  of  intervals,  0  =  no  overlap).  Then,  given 
repeated  application  of  our  candidate  procedures,  such  as  we  have  done  in  Section  5.5, 
the  proportion  of  time  that  Overlap  =  1  occurred  is  an  estimate  of  P*. 

5.3.  Candidate  Results:  One  Example 

In  this  section,  we  compare  the  results  from  one  example  of  using  Candidate  SI 

and  Candidate  S2,  using  data  from  Scenario  1  (AAAA).  When  a  candidate  determines 

that  the  stopping  criterion  has  been  met,  we  say  that  the  candidate  issues  a  stop  order. 

Figure  5.1  presents  a  visual  display  of  the  results  of  using  SI,  while  Table  5.2  gives  a 

summary  of  Si’s  performance.  Likewise,  Figure  5.2  presents  a  similar  visual  display  for 

a  single  realization  of  using  S2  on  the  same  data  and  Table  5.3  presents  a  summary  of 

S2’s  performance.  Since  S2  involves  the  stochastic  bootstrap  sampling  procedure,  this 
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realization  may  vary  each  time  the  procedure  is  performed.  The  particular  realization  of 
S2  that  appears  in  Figure  5.2  and  Table  5.3  was  chosen  for  display  because  in  this 
realization,  N2  ~E(N2).  Specifically,  N2  =  97  in  the  chosen  realization,  and  our 
estimate  of  E (N2 )  is  99.4.  The  variability  of  S2  when  performed  on  Scenarios  1-3  is 
assessed  in  Section  5.4. 

In  both  Figure  5. 1  and  Figure  5.2,  when  both  confidence  interval  bands  cross  the 
relative  error  threshold  bands,  the  candidate  issues  a  stop  order.  From  these  figures,  we 
observe  that  there  are  two  subtle  differences  in  these  processes,  which  are  noted  below: 

1.  The  process  used  by  both  candidates  requires  that  both  the  upper  and  lower 
confidence  bands  be  within  the  relative  error  threshold  bands  before 
determining  that  enough  runs  have  been  made.  In  Candidate  Si’s  case,  the 
confidence  intervals  are  symmetric;  therefore,  on  any  given  run,  either  neither 
or  both  of  the  confidence  bands  are  within  the  relative  error  threshold  bands. 
This  is  not  the  case  with  Candidate  S2.  Since  bootstrap  percentile  intervals 
are  not  necessarily  symmetric,  on  any  given  run,  zero,  one,  or  two  confidence 
bands  may  be  within  the  relative  error  threshold  bands. 

2.  If  we  allow  the  BRAWLER  runs  to  continue  after  a  stop  order  is  issued,  we 
see  that  Candidate  SI  is  more  consistent  in  issuing  another  stop  order  on  the 
next  run  after  the  first  stop  order  was  issued.  S 1  will  only  reverse  its  decision 
to  issue  a  stop  order  if  the  variance  of  the  sample  data  increased  sufficiently 
(due  to  the  new  data)  from  the  value  calculated  on  the  previous  run.  In 
contrast,  in  a  similar  situation  it  is  quite  possible  that  S2  will  reverse  its 
decision  to  issue  a  stop  order  on  the  next  run  even  if  the  variance  of  the 
sample  data  did  not  increase.  This  phenomenon  is  due  to  the  stochastic  nature 
of  the  bootstrap  sampling  process. 

Figure  5.1  and  Figure  5.2  allow  us  to  compare,  for  each  run,  the  lower  and  upper 
thresholds  that  meet  our  stopping  criterion  to  the  lower  and  upper  endpoints  of  our 
confidence  interval. 
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Figure  5.1.  Candidate  SI,  Scenario  1 


Table  5.2.  Summary  of  Candidate  SI  Performance,  Scenario  1 


#  Runs  to  Stop  (N\) 

122 

Estimated  Mean  (E(AKPFi)) 

32.39 

90%  Standard  Cl 

(29.45,  35.32) 

Cl  Length 

5.87 

Overlap 

1 

Using  Candidate  SI,  the  stopping  criterion  was  met  after  122  runs  were 
performed.  At  122  runs,  the  point  estimator  for  the  mean  was  32.39  and  the  associated 
90%  standard  confidence  interval  was  (29.45, 35.32),  which  has  a  length  of  5.87.  The 
standard  confidence  interval  from  122  runs  overlaps  the  Truth  model ,  the  90% 
confidence  interval  at  200  runs,  which  is  (30.07,  34.28). 
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Figure  5.2.  Candidate  S2,  Scenario  1 


Table  5.3.  Summary  of  Candidate  S2  Performance,  Scenario  1 


#  Runs  to  Stop  (N2) 

97 

Estimated  Mean  (E(AKPF2)) 

33.92 

90%  Standard  Cl 

(30.67,  37.13) 

Cl  Length 

6.47 

Overlap 

1 

Using  Candidate  S2,  the  stopping  criteria  for  the  particular  realization  shown  in 
Figure  5.2  was  met  after  97  runs  had  been  performed.  At  97  runs,  the  point  estimator  for 
the  mean  was  33.92  and  the  associated  90%  bootstrap  percentile  interval  was 
(30.67,  37.13),  which  has  a  length  of  6.47.  The  bootstrap  percentile  interval  from  97  runs 
also  overlaps  the  Truth  model. 
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5.4.  Variability  of  Candidate  S2 


We  have  seen  that  the  bootstrap  method  itself  is  a  stochastic  simulation  process. 
Because  of  this,  each  time  we  use  Candidate  S2,  we  may  obtain  a  different  answer  for 
both  the  number  of  runs  required  to  meet  the  stopping  criterion  and  the  value  of  AKPF 
that  is  returned.  Table  5.4  displays  the  results  from  20  replications  of  S2  on  each  of  our 
three  scenarios.  To  aid  making  a  comparison.  Table  5.5  shows  the  similar  results  from 
using  Candidate  SI. 


Table  5.4.  Results  from  20  Replications  of  Candidate  S2 


Candidate  S2 

Scenario  3 

Replication 

n2 

AKPF 

n2 

AKPF 

1 

103 

33.78 

100 

26.93 

123 

28.80 

2 

98 

33.95 

88 

27.81 

120 

28.22 

3 

100 

33.80 

86 

27.93 

137 

29.56 

4 

106 

33.70 

97 

27.36 

131 

29.61 

5 

99 

34.04 

96 

27.18 

118 

27.65 

6 

99 

34.14 

86 

27.97 

139 

30.19 

7 

102 

33.84 

93 

27.44 

134 

29.48 

8 

100 

33.84 

99 

26.92 

127 

29.25 

9 

102 

33.63 

27.32 

136 

29.15 

95 

33.69 

26.96 

129 

28.92 

11 

100 

33.69 

84 

27.28 

125 

28.94 

12 

83 

34.40 

84 

27.25 

123 

28.79 

13 

106 

33.70 

125 

29.19 

14 

102 

33.51 

93 

27.49 

129 

29.08 

15 

102 

33.73 

84 

27.29 

126 

29.60 

16 

96 

33.14 

92 

27.36 

123 

28.66 

95 

33.80 

99 

26.84 

121 

28.60 

18 

99 

33.99 

96 

131 

29.50 

19 

99 

34.15 

mm 

27.34 

127 

29.33 

102 

33.46 

26.88 

125 

28.86 

Expected 

Value 

99.4 

33.80 

127.45 

29.07 
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Table  5.5.  Candidate  SI  Results 


Candidate  SI 

Scenario  1 

Scenario  2 

Scenario  3 

N\  AKPF 

V,  AKPF 

V,  AKPF 

Value 

122  32.39 

117  27.85 

190  28.53 

We  see  that  on  average,  Candidate  S2  required  18.5%  less  runs  than  Candidate  SI  to 
issue  a  stopping  order  on  Scenario  1,  20.1%  less  on  Scenario  2,  and  32.9%  on 
Scenario  3.  These  results  are  displayed  in  Table  5.6. 


Table  5.6.  Results  of  SI  and  S2  on  Scenarios  1-3 


Candidate  Si 

Scenario  1 

Scenario  2 

Scenario  3 

mmm 

E(Ni)  E(AKPFi) 

122.0 

32.39 

117.0  27.85 

190.0  28.53 

WBM 

99.4 

33.80 

93.5  27.27 

127.5  29.07 

%  Runs  Saved 
with  S2 

18.5 

20.1 

32.9 

5.5.  Effect  of  the  Order  of  BRAWLER  Data 
5.5.1.  Introduction 

Since  both  candidates  are  sequential  procedures,  the  order  in  which  the 
BRAWLER  data  is  obtained  affects  both  when  a  stop  order  is  issued  and  the  estimate  of 
expected  AKPF  that  is  returned.  When  we  mention  the  order  in  which  the  data  was 
obtained,  we  are  referring  to  the  replication  number  on  which  a  specific  (multivariate) 
BRAWLER  data  is  obtained.  In  our  case,  we  are  interested  in  the  bivariate  data 
consisting  of  APKF  and  Firings  values.  If  for  some  reason  the  data  obtained  before  a  stop 
order  was  issued  is  biased,  our  sequential  candidate  procedures  return  biased  estimates  of 
the  expected  AKPF.  Therefore,  we  determined  it  is  necessary  to  investigate  the 
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performance  of  each  candidate  on  different  random  orderings  of  the  AKPF  data  obtained 
from  our  200  BRAWLER  runs. 

5.5.2.  Results  of  Investigation 

We  performed  both  candidate  procedures  on  20  different  randomized  orderings  of 
the  data  from  200  BRAWLER  runs  for  Scenario  1.  Each  time  we  perform  a  procedure  on 
a  different  ordering  of  the  data  we  say  that  we  are  performing  a  metareplication  of  the 
procedure.  Although  we  performed  20  replications  of  Candidate  S2  in  the  experiment 
mentioned  in  Section  5.4,  in  this  experiment  we  performed  only  a  single  metareplication 
of  S2  for  each  of  the  20  orderings  of  the  data. 

Table  5.6  summarizes  the  performance  of  each  of  our  candidates  on  the  different 
metareplications.  In  its  top  section,  the  table  shows  the  90%  confidence  interval  that  was 
constructed  for  our  Truth  model  for  comparison  purposes.  In  the  lower  section,  the  table 
lists  data  on  the  performance  of  each  candidate.  The  table’s  first  column  lists  the  labels 
for  the  different  orderings  of  BRAWLER  data  that  were  used.  These  include  the  original 
order  of  the  data  and  the  20  reorderings  (denoted  Reorder  r,  where  r  e  [1,20] ).  In  the 

successive  columns,  the  following  data  is  displayed  for  Candidate  i,  first  for  i  =  1  and 
then  for  i  =  2: 

1.  The  number  of  runs  M  (denoted  in  the  table  as  #Runs)  required  before 
Candidate  i  issued  a  stop  order, 

2.  The  AKPF  value  returned  by  Candidate  i  after  N\  runs, 

3.  The  associated  90%  confidence  interval  (Cl)  for  the  AKPF  in  number  2  above, 
which,  in  the  case  of  Candidate  S2  is  the  90%  bootstrap  percentile  interval 

(PI). 

4.  The  value  of  Overlap,  as  defined  in  Section  5.2. 
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Table  5.6.  Candidate  Performance  on  Different  Orderings  of  BRAWLER  Data 


#Runs 

AKPF 

90%  Cl 
for  AKPF 

Overlap 

#Runs 

AKPF 

90%  PI 
for  AKPF 

Overlap 

Truth  Model 

- 

- 

- 

200 

32.18 

30.071 

34.284 

- 

- 

- 

Candle 

late  1 

I  Candidate  2  -  Single  Realizai 

tion 

Original  Order 

122 

32.39 

29.45 

35.32 

1 

30.667 

37.13 

1 

Reorder  1 

119 

30.86 

28.07 

33.64 

1 

31.27 

28.191 

34.36 

1 

Reorder  2 

120 

31.07 

28.26 

33.88 

1 

85 

32.76 

29.133 

35.95 

1 

Reorder  3 

117 

31.49 

28.63 

34.35 

1 

87 

30.61 

27.886 

33.36 

1 

Reorder  4 

120 

33.85 

30.79 

36.90 

1 

94 

34.71 

31.402 

37.79 

1 

Reorder  5 

115 

33.81 

30.75 

36.87 

1 

99 

33.67 

30.527 

36.92 

1 

Reorder  6 

119 

30.14 

27.41 

32.86 

1 

102 

28.89 

26.046 

31.75 

1 

Reorder  7 

140 

30.21 

27.47 

32.96 

1 

94 

29.54 

1 

Reorder  8 

106 

35.82 

32.58 

39.06 

1 

91 

36.17 

1 

Reorder  9 

86 

33.46 

30.43 

36.48 

1 

75 

32.07 

29.091 

35.24 

1 

Reorder  10 

119 

29.58 

26.90 

32.25 

1 

62 

33.82 

30.476 

36.94 

1 

Reorder  1 1 

138 

33.08 

30.09 

36.07 

1 

122 

33.31 

30.023 

36.60 

1 

Reorder  12 

98 

33.02 

30.02 

36.02 

1 

34 

42.00 

37.865 

46.11 

0 

Reorder  13 

105 

32.60 

29.67 

35.53 

1 

83 

32.30 

29.742 

35.37 

1 

Reorder  14 

141 

30.96 

28.15 

33.77 

1 

111 

30.59 

28.004 

33.46 

1 

Reorder  15 

122 

33.32 

30.30 

36.34 

1 

95 

32.76 

29.552 

35.96 

1 

Reorder  16 

129 

32.29 

29.37 

35.22 

1 

97 

30.85 

27.83 

33.55 

1 

Reorder  17 

100 

32.42 

29.53 

35.30 

1 

78 

31.35 

28.417 

34.45 

1 

Reorder  1 8 

101 

30.33 

27.62 

33.05 

1 

81 

28.45 

25.603 

31.28 

1 

Reorder  19 

106 

32.28 

29.37 

35.19 

1 

84 

31.41 

28.464 

34.46 

1 

Reorder  20 

121 

31.55 

28.73 

34.37 

1 

106 

30.15 

27.145 

33.16 

1 

Mean  Value 
Overlap  Prop 

116.38 

32.12 

- 

1.00 

88.10 

32.41 

- 

- 

0.95 

5.5.3.  Expected  Number  of  Runs 

We  constructed  approximate  90%  standard  confidence  intervals  for  the  estimated 
means  of  the  number  of  runs  N\  and  N2  required  before  a  stop  order  was  issued.  These 
intervals  appear  in  Figure  5.1below.  From  this  figure  we  may  infer  that  E(A^)  <  E(M); 
however,  a  formal  comparison  is  made  in  Section  5.5.5.  Our  estimate  of  E(M)  based  on 
the  sample  data  is  1 16.38,  which  means  that  we  would  expect  to  perform  approximately 
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117  runs  of  Scenario  1  before  meeting  our  stopping  criterion  using  Candidate  SI. 
Furthermore,  we  state  with  90%  confidence  that  E(N\)  lies  in  the  interval  (111  .04, 
121.73).  In  contrast,  our  estimate  ofE(N2)  is  88.10,  so  we  would  expect  to  perform 
approximately  89  runs  of  Scenario  1  before  stopping  when  using  Candidate  S2,  and  we 
state  with  90%  confidence  that  E(N2)  lies  in  the  interval  (81.15,  95.04).  So,  based  on  20 
metareplications  of  Scenario  1,  we  see  that  on  average  S2  issued  a  stopping  order  using 
24%  fewer  runs  than  SI. 


130  -I 

120  - 

j  121.73 

[  W  116.38 

«  HO- 

A  111.04 

(E  100  - 

- 

=tfc 

T  95.04 

90  - 

■  88.10 

80- 

70  - 

A  81.15 

- 1 - 

1  2 
Candidate 


Figure  5.1.  90%  Confidence  Intervals  for  E(N,)  and  E(N7) 

5.5.4.  Expected  Value  of  AKPF 

Next,  we  constructed  approximate  90%  standard  confidence  intervals  for  the 
estimated  mean  of  AKPF\  (the  AKPF  resulting  from  a  stop  order  issued  by  SI  at  run  N\), 
and  AKPFi  (the  AKPF  for  S2  at  run  N2).  These  intervals  appear  in  Figure  5.2.  From  this 
figure  we  may  infer  that  E(AKPF\)  =  E(AKPF2)\  however,  the  formal  comparison  appears 
in  Section  5.5.5.  Our  point  estimate  of  E(AKPF\ )  based  on  the  sample  data  is  32. 12. 
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sample  data  is  32.12.  Furthermore,  we  state  with  90%  confidence  that  E(AKPFi)  for 
Scenario  1  (when  our  stopping  criterion  is  met  using  SI)  lies  in  the  interval  (31.54, 
32.70).  Our  point  estimate  of  E(AKPF2)  is  32.41.  we  can  state  with  90%  confidence  that 
E{AKPF2)  for  Scenario  1  (when  our  stopping  criterion  is  met  using  S2)  lies  in  the  interval 
(31.30, 33.51).  The  confidence  interval  for  S2  is  wider  than  that  for  SI  due  to  the 
variability  introduced  by  the  stochastic  nature  of  the  bootstrap,  as  discussed  in  Section 
5.4.  For  comparison  purposes.  Figure  5.3  displays  the  90%  confidence  interval 
constructed  from  our  Truth  model .  This  figure  appears  here  to  show  that  both  point 
estimates  for  E(AKPF\)  and  E(AKPF2 )  compare  well  with  the  point  estimate  of  the 
expected  AKPF  from  our  Truth  model .  However,  the  reader  should  not  directly 
compare  the  length  of  the  intervals  in  Figure  5.2  and  Figure  5.3.  The  interval  in  Figure 
5.3  is  used  to  estimate  the  expected  AKPF  of  Scenario  1  based  on  a  sample  size  of  200 
BRAWLER  replications.  In  contrast,the  intervals  in  Figure  5.2  are  used  to  estimate  the 
mean  AKPF  response  for  Scenario  1  that  is  obtained  by  each  candidate  based  on  a  sample 
size  of  20  metareplications  of  the  200  BRAWLER  replications. 
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Figure  5.2,  90%  Confidence  Intervals  for  E (AKPF,)  and  E (AKPF,) 


Note:  We  cannot 
directly  compare  the 
length  of  the  intervals 
in  Figure  5.2  and 
Figure  5.3.  See  text 
above  for  an 
explanation. 


T ruth  Model 


Figure  5.3. 90%  Confidence  Interval  for  E (AKPF)  from  Truth  Model 


5.5.5.  Paired  Comparisons 

We  have  been  comparing  the  two  candidates  on  the  basis  of  two  different 
measures,  namely: 
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1.  The  number  of  runs  required  before  a  stop  order  is  issued  (a  MOM),  and 

2.  The  value  of  AKPF  returned  and  its  associated  confidence  interval. 

We  may  make  these  comparisons  by  constructing  a  90%  confidence  interval  for  the 
difference  in  the  expected  values  in  numbers  1  and  2  above  for  each  candidate,  which  we 
may  write  explicitly  as: 

3.  E(Nt)-E(N2) 

4.  E(  AKPF, )  -  E(  AKPF2 ) 

Since  we  can  pair  the  observations  in  Table  5.6  according  to  the  ordering  of  the  data  they 
were  subjected  to,  we  may  make  a  paired  comparison  to  estimate  the  values  in  numbered 
items  3  and  4  above.  We  may  do  this  by  constructing  a  paired-t  confidence  interval  for 
the  difference  in  the  expected  values  we  are  interested  in.  A  90%  paired-?  confidence 
interval  for  E ( A, )  -  E(  N2 )  appears  in  Figure  5.4.  In  this  figure  we  can  see  that  the 
difference  in  expected  values  of  N,  and  N2  is  significant  at  a  =  .10  since  the  confidence 

interval  does  not  include  zero.  The  interval  displayed  is  (22.91,  33.66)  and  the  point 
estimator  of  E(  N] )  -  E(  N2 )  is  28.29. 

Likewise,  a  paired-t  90%  confidence  interval  for  E(  AKPF, )  -  E(  AKPF2  )appears 
in  Figure  5.5.  In  contrast  to  Figure  5.4,  Figure  5.5  shows  that  the  difference  in  expected 
values  of  AKPFX  and  AKPF2  is  not  significant  at  a  =  .10  since  the  portrayed  confidence 

interval  includes  zero.  The  interval  displayed  is  (-1.20, 0.62)  and  the  point  estimator  of 
E(  AKPF, )  -  E(  AKPF2 )  is  -0.29. 
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Figure  5.4.  90%  Paired-t  Confidence  Interval  for  E(/V,)  -  E(N2) 


Figure  5.5.  90%  Paired-t  Confidence  Interval  for  E (AKPF,)  -  E (AKPFJ 


5.5.6.  Probability  of  Qverlai 


As  can  be  seen  in  Table  5.6,  the  proportion  of  time  the  confidence  intervals 
constructed  using  Candidate  i  and  the  confidence  interval  constructed  from  our  Truth 
model  overlapped  was  1.00  for  i  =  1  and  0.95  for  i  =  2.  We  use  these  values  to  estimate 
Pi  and  P2,  respectively.  So,  based  on  our  sample  of  20  metareplications,  we  estimate  the 
value  of  Pi  to  be  1 .00  and  P2  to  be  0.95.  We  constructed  a  90%  confidence  interval  for 
E(Pi  -  P2)  which  appears  in  Figure  5.6.  Figure  5.6  shows  that  the  expected  value  of  Pi  - 
P2  is  not  significant  at  a  -  .10  since  the  90%  confidence  interval  includes  zero.  The 
interval  displayed  is  (-0.03,  0.13)  and  the  point  estimator  of  Px  -  P2  is  0.05. 


Figure  5.6.  90%  Confidence  Interval  for  E(Pi  -  P2) 


6.  Conclusions 


6.1.  Summary 

We  have  demonstrated  proof  of  the  SIMAF  concept  through  our  pilot  program 
that  demonstrates  the  cross-resolution  modeling  process  through  the  calibration  of  the 
campaign-level  model  THUNDER  with  the  engagement-level  model  BRAWLER.  Since 
both  of  these  models  are  recommended  members  of  the  USAF  Analysis  Toolkit,  they  will 
continue  to  be  used  for  several  years  while  the  USAF  transitions  to  newer  models.  Using 
cross-resolution  modeling,  we  can  ensure  that  outputs  from  THUNDER  are  consistent 
with  their  BRAWLER  counterparts.  We  can  also  ensure  that  a  pedigree  is  established  for 
data  crossing  levels  of  model  resolution.  We  have  shown  how  a  particular  MOE  from 
BRAWLER,  the  adjusted  kills  per  firing  (AKPF),  is  transformed  into  input  data  for 
THUNDER  air-to-air  probability  of  kill  (PK).  The  procedures  in  Chapter  5  are  different 
ways  of  estimating  the  expected  value  of  the  AKPF  in  BRAWLER,  which  is  the  value 
that  is  used  as  the  THUNDER  air-to-air  PK.  All  of  the  procedures  presented  in  Chapter  5 
provided  sufficiently  accurate  estimates  of  E(AKPF). 

Using  sequential  procedures,  we  can  minimize  the  number  of  BRAWLER  runs 
necessary  to  provide  us  with  a  sufficiently  accurate  estimate  of  E {AKPF).  We  define 
sufficient  accuracy  up  front  in  terms  of  a  relative  error  e  and  significance  level  a,  which 
determines  our  stopping  criterion.  A  (1  -a)%  confidence  interval  may  be  automatically 
constructed  after  each  run  that  lets  us  determine  if  we  have  met  our  stopping  criterion. 
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Minimizing  the  number  of  runs  in  this  way  allows  us  to  increase  the  speed  of  the 
THUNDER  air-to-air  calibration  process. 

In  particular,  we  found  that  the  sequential  Candidate  SI,  which  uses  the  bootstrap, 
performs  better  than  the  sequential  Candidate  S2,  which  uses  standard  normal  theory.  In 
our  20  metareplications  of  Scenario  1,  we  found  that  on  average  S2  issued  a  stopping 
order  using  24%  fewer  runs  than  SI,  while  returning  statistically  similar  results  for  both 
the  expected  value  of  AKPF  and  the  probability  of  overlapping  the  confidence  interval 
from  our  Truth  model .  However,  both  procedures  perform  well,  so  they  are  both 
recommended  for  SIMAF  use. 

6.2.  Discussion 

It  should  be  noted  that  bootstrap  methods  are  an  alternative  to  standard  statistical 
procedures,  not  a  replacement  for  them  (Cheng,  1995).  Bootstrap  methods  used  in 
conjunction  with  classical  statistical  procedures  provide  both  the  analyst  and  the  decision¬ 
maker  a  great  deal  of  insight  about  the  distribution  of  a  statistic  of  interest. 

Through  the  use  of  bootstrap  sampling,  we  found  that  a  major  asset  of  the 
bootstrap  is  its  ability  to  build  up  a  picture  of  the  distribution  of  a  statistic  of  interest 
based  on  a  sample.  This  is  especially  valuable  for  smaller  sample  sizes  where  a  normality 
assumption  may  not  be  appropriate.  Not  only  is  the  bootstrap  a  powerful  tool,  but  it  can 
be  implemented  easily  and  quickly  by  performing  a  fairly  trivial  computer  simulation. 

The  simulation  can  be  written  in  any  number  of  programming  languages  -  general 
programming  languages  such  as  C++  or  FORTRAN,  or  mathematical  programming 
languages  are  examples. 
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With  modem  mathematical  programming  languages,  such  as  Mathcad,  the  process 
of  building  a  picture  of  the  distribution  of  a  statistic  can  also  be  graphically  animated. 
Using  the  Mathcad  6.0  program  we  created,  we  were  able  to  do  this  quite  effectively. 
However,  it  took  a  longer  time  to  learn  how  to  exploit  this  capability  than  it  did  to 
implement  the  actual  bootstrap  procedure.  Yet  this  capability  is  especially  valuable  in 
this  day  and  age  where  emphasis  on  the  visual  display  of  data  is  ever  increasing.  Data 
visualization  techniques  are  a  key  part  of  the  SIMAF  concept,  therefore  the  bootstrap  is  a 
very  appropriate  tool  for  SIMAF  use. 

6.3.  Additional  Research 

Another  useful  tool  for  the  simulation  analyst  is  a  between-run  prediction  of  the 
number  of  runs  that  must  still  be  performed  (runs  to  go)  in  order  to  meet  the  chosen 
stopping  criterion.  Some  initial  experimentation  in  this  area  was  accomplished  in  this 
research,  but  due  to  time  limitations  was  not  fully  explored;  therefore,  only  initial  results 
are  presented  here.  We  considered  two  candidates  PI  and  P2  for  predicting  the  number  of 
runs  to  go.  These  candidates  are  closely  related  to  candidates  SI  and  S2  for  determining 
if  our  stopping  criterion  was  met. 

6.3.1.  Candidate  PI 

Candidate  PI  makes  use  of  another  procedure  presented  by  Law  and  Kelton 
(1991)  in  order  to  make  such  a  prediction  using  the  standard  confidence  interval.  We 
combined  this  procedure  with  Law  and  Kelton’ s  sequential  procedure  described  in 
Section  4.6.2.  Using  this  combination,  we  were  able  to  obtain  a  prediction  of  the  number 
of  runs  yet  to  be  performed  after  every  run  that  did  not  meet  the  stopping  criterion. 
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Law  and  Kelton  present  an  approximate  expression  for  the  number  of  replications 
needed  to  obtain  a  relative  error  y ,  given  by: 


»*(y)  =  nmi -I 


>  n .  hthtZlR '  <  y'| 

|x(n)| 


(6.1) 


where  t(i_u_a/2)  •  ■y]d\n)/i  =  hl(i,a)  is  the  same  confidence  interval  half-length  as  in 


(4.21)  except  that  it  is  now  based  on  i  >  n  observations  instead  of  n.  We  make  an 
assumption  that  our  estimates  of  the  population  mean  and  population  variance, 

A  /s,  0 

x(n )  =  6{n)  and  a  (n) ,  respectively,  do  not  change  as  more  replications  are  performed. 
Then  n*(y)  is  approximated  by  the  smallest  integer  i  that  satisfies 
i  >  <72(n)[z(]_a/2)/y'  •  x(n)f .  Using  (6.1),  we  were  able  to  obtain  a  prediction  of  the 

number  of  runs  yet  to  be  performed  after  every  run  that  did  not  meet  the  stopping 
criterion.  We  decided  that  the  best  way  to  judge  the  performance  of  this  procedure  was  to 
compare  the  following: 

1.  A  plot  of  the  predicted  number  of  runs  to  go  versus  the  current  run  number. 

2.  A  plot  of  the  actual  number  of  runs  to  go  versus  the  current  run  number.  This 
is  a  straight  line  we  denote  as  the  Truth  Line. 

Based  on  a  visual  inspection  of  these  plots,  we  conclude  that  this  procedure  perforins  well 
on  our  three  scenarios.  More  formal  methods  forjudging  performance,  such  as 
calculating  a  mean  squared  error  based  on  the  deviations  from  the  Truth  Line  could  be 
implemented  if  so  desired.  Our  results  are  presented  in  Figure  6.1  through  Figure  6.3. 
Since  we  use  a  set  of  20  pilot  runs  to  obtain  an  initial  prediction,  predictions  for  runs 
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Figure  6.1.  Predictions  of  Runs  to  Go  for  Scenario  1  using  PI 


Figure  6.2.  Predictions  of  Runs  to  Go  for  Scenario  2  using  PI 
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Figure  6.3.  Predictions  of  Runs  to  Go  for  Scenario  3  using  PI 

1-19  are  not  made,  so  for  each  of  these  runs  the  predicted  number  of  runs  to  go  is 
assigned  a  value  of  zero. 

6.3.2.  Candidate  P2 

The  bootstrap  procedure  is  normally  involved  with  generating  bootstrap  samples 
and  then  obtaining  bootstrap  statistics  for  each  of  these  bootstrap  samples.  A  key  aspect 
of  the  bootstrap  process  is  this:  Given  a  statistic  of  interest  for  a  sample,  such  as  the 
mean,  a  large  sample  of  B  bootstrap  statistics  has  nearly  the  same  variance  as  the  original 
statistic  of  interest.  Therefore,  the  variance  of  the  bootstrap  statistics  can  be  used  to 
estimate  the  variance  of  the  original  statistic  (Cheng,  1995). 
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We  estimate  a  standard  error  (standard  deviation  of  the  sample  mean)  with  the 
standard  deviation  of  a  large  sample  of  B  bootstrap  means,  denoted  6§.  ( B ) .  We  write 

this  as  -y/tr 2(n)/n  *  6 (b).  The  bootstrap  estimate  of  standard  error  is 

l 

<v  ( B )  =  jx  (b) -e*}j(B- 1)  j 2  (6.2) 

where  6 * (b)  is  the  estimate  of  6  for  the  bootstrap  sample  x*b ,  <7*-  is  the  estimated 

standard  error  of  6*  for  the  bootstrap  sample  jc**  and  0*=  (l/fi)  ■  ^fh_^*(b)  . 

We  use  the  bootstrap  estimate  of  standard  error  in  place  of  the  original  estimate  of 
standard  error  to  create  a  slight  variation  on  (6.1),  as  follows: 


nr(y)  =  rmn< 


i>n: 

^(i-l.l-a/2)  '  (#)■  Vrt/Vi) 

r\rf  > 

;  |*00| 

-7 

(6.3) 


We  then  sequentially  predict  the  number  of  runs  to  go  after  each  run  is  completed 
using  (6.3).  The  reason  for  using  (6.3)  instead  of  (6.1)  is  to  investigate  if  using  the 
bootstrap  estimate  of  standard  error  yields  different  prediction  results  than  using  the 
standard  estimate  does.  The  results  for  P2  used  on  Scenario  1  are  presented  in  Figure  6.4 
for  B  =  501.  The  series  denoted  Rep  i,  where  i  e  [1,10] ,  plots  the  predictions  of  the 
number  of  runs  to  go  on  the  y-axis  against  the  number  of  runs  completed  on  the  x-axis  for 
each  of  10  replications  of  P2.  The  series  denoted  Mean  Predict  to  Go  plots  the  mean  of 
the  predictions  made  after  each  run  is  completed  for  10  replications  of  P2.  Finally,  the 
series  denoted  Mean  Actual  to  Go  plots  the  actual  number  of  runs  to  go,  based  on  using 
run  121  as  the  final  run  (stopping  run).  Run  121  was  chosen  because  for  the  10 


6-7 


replications  of  P2,  the  mean  number  of  runs  to  meet  the  stopping  criterion  is  121.2.  Table 
6.1  presents  the  results  of  10  replications  of  P2  on  Scenario  1,  using  the  convention  in 
Table  5.4. 

Initial  results  seem  to  indicate  that  there  is  a  difference  between  the  predictions 
made  by  PI  and  the  mean  of  the  predictions  made  by  P2  on  a  given  run.  It  appears  that 
the  mean  of  predictions  made  by  P2  may  perhaps  be  slightly  more  accurate  than  a 
prediction  made  by  PI  on  runs  distant  from  the  stopping  run.  However,  the  predictions 
made  by  P2  on  a  given  run  have  a  very  high  variance,  most  noticeably  when  the  given  run 
is  close  to  the  stopping  run.  Again,  this  is  due  to  the  stochastic  nature  of  the  bootstrap. 

Therefore,  based  on  our  initial  experimentation,  a  single  prediction  of  the  number 
of  runs  to  go  made  by  P2  using  the  bootstrap  estimate  of  standard  error  does  not  appear  to 
be  useful.  However,  if  the  prediction  process  is  replicated,  the  mean  of  the  predictions 
made  by  P2  appears  to  be  fairly  accurate  on  runs  far  from  the  stopping  run.  Since  P2’s 
prediction  process  may  be  replicated  quickly,  P2  may  be  a  useful  tool  for  predicting  the 
number  of  runs  to  go.  The  best  strategy  may  be  to  rely  on  P2  for  predictions  when  the 
current  run  is  far  from  the  stopping  run,  and  then  start  relying  on  PI  when  the  predictions 
from  P2  become  sufficiently  small  (for  example,  30  runs  to  go). 
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Appendix  A.  bootstrav.mcd 


BOOTSTRAP  SAMPLING  PROGRAM  using  Mathcad  6.0 

(Based  on  RAND  Corp.'s  enoughruns.F  Fortran  program  designed 
to  determine  when  enough  BRAWLER  runs  have  been  made) 

GLOBAL  PARAMETERS: 

B  :  —501  B  is  the  #  of  bootstrap  reps.  WARNING:  INCREASING  B  INCREASES 

B  can  be  any  number  you  PROCESSING  TIME  EXPONENTIALLY 
choose.  Typically,  B  is 
chosen  such  that 
500<B<2000. 

In  enoughruns.F ,  B:=501 

p  and  q  are  used  to  look  at  the  advantages  of  using  bootstrapping  with 
less  runs  vs.  normal  methods  with  more  runs. 

p  :-94  p  is  the  number  of  data  points  in  the  output  data  vector  (called  Data 
below)  that  you  would  like  to  consider  for  bootstrapping  (1<p<q<=k, 
where  k  is  the  total  number  of  data  points  you  have  available). 

Typically,  you  would  never  set  p<20  when  using  bootstrapping. 

In  enoughruns.F,  p>19,  since  20  is  the  minimum  number  of  BRAWLER 
runs  that  must  be  made  before  enoughruns.F  begins  bootstrapping  runs. 

q  :=1 17  q  is  the  number  of  data  points  in  the  output  data  vector  (called  Data 

below)  that  you  would  like  to  consider  for  non-bootstrapping  (1<p<q<=k, 
where  k  is  the  total  number  of  data  points  you  have  available). 

a  ’ 1 0  a  is  the  significance  level  (0<a<1 ) 

Y  :=.10  y  is  the  relative  error  (0<y<1) 

When  we  obtain  an  estimate  of  p,  the  population  mean,  with  a  relative 
error  of  y  and  a  confidence  level  of  1 00(1  -a)%,  we  say  that  we  have 
enough  simulation  data,  or,  in  other  words,  "enough  runs." 

Notes:  ORIGIN  has  been  set  to  1  in  Menu  Item  Math  -  Built-In  Variables 

Also  on  Menu:  Math  -  Automatic  Mode  must  be  enabled  for  animation  to  work 
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Function  sumfirings  counts  the  total  # 
of  firings  in  all  the  runs  in  X.  sumfirings 
expectsX  to  be  a  matrix  of  AKPFs  and 
Firings. 


Function  weightPks  creates  an  equally  weighted 
vector  (w.r.t.  Firings)  of  AKPFs  suitable  for  boot¬ 
strapping.  The  AKPF  for  a  run  is  entered  in  the 
vector  a  number  of  times  =  to  the  #  of  Firings  in 
that  run 


sumfirings(X) 


sum<-0 
r«-rows(X) 
for  i  e  1 ..  r 
sumf-sum-f-Xj  2 
sum 


weightPks(X)  !- 


s*-sumfirings(X) 

r«-rows(X) 

Wj<-0 


for  i  e  1 ..  r 
if  X;  2>0 

for  j  e  1 ..  Xj  2 

Wt^Xi, 

W«-stack(W,  Wt) 

Wt«-0 

W*-  submatrix  (w,  2 ,  s  -f- 1 , 1 , 1 ) 
W 


Function  bootstraps  bootstraps  1  value 
from  the  runs  data  (original  sample).  To 
obtain  a  bootstrap  sample,  the  function 
must  be  replicated  k  times  where  k  is  the 
size  of  the  original  sample.  It  takes  a 
vector  X  of  size  k  and  returns  1  randomly 
sampled  element  of  X 


Function  bootstraps  pair  bootstraps  1  ordered 
pair  (AKPF,  Firings)  from  the  runs  data  (original 
sample).  To  obtain  a  bootstrap  sample,  the  function 
must  be  replicated  N  times  where  N  is  the 
size  of  the  original  sample.  It  takes  a  matrix  X  of  size 
Nx2  and  returns  1  randomly  sampled  ordered  pair 
(AKPF, Firings) 


boots  trap_l(X) 


r«-rows(X) 
N«— X 
uf-md(r) 
u«-ceil(u) 

m_Nu 

x 


bootstrap_lpair(  X) 


r«-rows(X) 

N*-X 

u«-md(r) 

u«-ceil(u) 

x.(N1)<U> 


Function  bootstrap_k  bootstraps  k 
values  from  an  original  sample  of  size  k. 

It  takes  a  vector  X  of  size  k  and  returns  a 
randomly  resampled  (bootstrap)version  of  X. 


Function  interval  takes  a  data  point 
x  and  a  relative  error  y  and  calculates 
the  interval  that  is  within  the  relative  error 
of  x. 


bootstrap Jc(  X) 


r<-rows(X) 

N*-X 


for  is  1  ..r 

|u«-md(r) 
uf-ceil(u) 


X 


interval  (x,y) 


interval  j  4- x-  (l  -  y) 
interval 2<-x*  (l  -f-y) 
interval 
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Function  stdev  calculates  standard  deviation. 

It  takes  a  vector  X  and  returns  the  standard 
deviation  of  X.  In  this  function,  the  variance  of 
X  is  modified  from  Mathcad's  built  in 
function  var,  since  Mathcad  uses  the  n 
divisor  instead  of  n-1. 


Function  stderr  calculates  standard  error. 

It  takes  a  vector  X  and  returns  the  standard 
error  of  X.  In  this  function,  the  variance  of 
X  is  modified  from  Mathcad's  built  in 
function  var,  since  Mathcad  uses  the  n 
divisor  instead  of  n-1 . 


stdev(X)  !- 


n«-rows(X) 

sd«-  Ivar(X) — — 
V  n-1 

sd 


stderr(X)  !- 


Function  bootstrap_mean  takes  a  vector  X, 
forms  one  bootstrap  sample  from  X  and 
returns  the  mean  of  the  bootstrap  sample 


Function  bootstrap_mean  takes  a  vector  X,  forms 
one  bootstrap  sample  from  X  and  returns  the  mean 
and  standard  error  of  the  bootstrap  sample 


bootmean(X) 


r«-rows(X) 
for  i  €  1 ..  r 
N* «-  bootstrap_  1  ( X ) 

m«-mean(  N) 
m 


bootstrap_mean_se(X)  I- 


r«-rows(X) 
for  i  g  1 ..  r 
hi*- boots  trap_l(X) 

m*-mean(N) 

sef-stderr(N) 


Function  Clhalflength  takes  a  vector  X  and 
a  significance  level  a  and  returns  the  half- 
length  of  a  100(1-a)%  confidence  interval 


Function  Clhalflength  takes  a  vector  X  and 
a  significance  level  a  and  returns  a  100(1-a)% 
confidence  interval  for  the  mean  of  X 


Clhalflength  (x,  a) 


n*-rows(X) 

se«-stderr(X) 

tvalue«-qt(  1  -  —  ,n-  1 

l  ^ 

Clhalflength*-  tvalue-se 
Clhalflength 


Cl(x,a) 


mn*-mean(X) 
hl«-  Clhalflength  (x ,  a) 
Clj^-mn-  hi 

Cl2«-mn-|-hl 

Cl 
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Function  percjnterval  takes  a  vector  X 
and  a  significance  level  a,  calculates 
the  5th  and  95th  percentile  of  X,  and 
returns  the  interval  having  those 
percentiles  as  endpoints. 


Function  Plhalflength  takes  a  vector  X 
and  a  significance  level  a  and  returns 
the  half-length  of  the  percentile  interval 
returned  by  function  percjnterval 


percjnterval (x, a)  •- 


r«-rows(X) 

Xs«-sort(X) 

Lo«-ceil  j  r— 

l  2, 


Upf-ceil 


r-  1. 


perc J  nterval  j 4-  Xs^Q 

perc Jnterva^^-  Xs^ 
percjnterval 


Plhalflength  (x ,  a) 


int«-  percjnterval  (x ,  a) 

inu  -  int , 

Plhalflength  < - - - - 

2 

Plhalflength 


Function  boot  _Cls  takes  a 
matrix  X  of  AKPF  &  Firings 
runs  data,  the  number  of 
bootstrap  replications  n, 
the  significance  level  a,  and 
a  kludge  factor,  kludge 
should  be  set  to  the  position 
of  the  first  nonzero  Firings 
entry  in  the  X  data  matrix. 

This  is  a  workaround  so  that 
if  initial  values  of  Firings  are 
equal  to  zero  the  program 
will  perform  correctly.  This 
function  returns  a  set  of  boot¬ 
strap  percentile  intervals  for  X. 
The  runs  data  is  read  in  one 
run  at  a  time  and  a  percentile 
interval  is  calculated  after  each 
new  run  is  added. 


boot_CIs  (x ,  n ,  a ,  kludge) 


q«-rows(X) 

Xs<-csort(x,  l) 
max«-Xs„  , 

qJ 

for  j  g  kludge.,  q 

Xj«-  submatrix  (x,  1  ,j ,  1 , 2) 
Xjw«-  weightPks(  Xj ) 
for  i  €  1 ..  n 
M-*-bootmean(  Xjw) 
Medj«-median(M) 

PI^  >  <-  percjnterval  (m  ,  a) 
Med 

<2> 


Function  boot  jneans  takes  a 
matrix  X  of  AKPF  &  Firings 
runs  data,  the  number  of 
bootstrap  replications  n, 
the  significance  level  a,  and 
a  kludge  factor,  kludge 
should  be  set  as  in  function 
boot_Cls.  This  function 
returns  a  set  of  n  bootstrap 
means  for  X,  where  the  boot¬ 
strap  mean  is  the  median  of 
the  means  of  each  bootstrap 
sample.  The  runs  data  is  read  in  one 
run  at  a  time  and  a  bootstrap 
mean  is  calculated  after  each 
new  run  is  added. 


boot_means 


(x,n, a, kludge)  '  - 


q«— rows(X) 
Xs«-csort(x,  l) 

max«— Xs  . 

q>i 


for  j  g  kludge.,  q 

Xj  submatrix  (x,  1  ,j,  1 ,2) 
Xjw«-  weightPks(  Xj ) 
for  i  g  1 ..  n 
M  j  4- bootmean(  Xj  w ) 


|  Medj  4—  median(  M ) 
Med 
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Function  boot_Z  takes  a  vector  X  and  the  number  of  bootstrap  replications  n  and  calculates  the 
bootstrap  version  of  the  Z  statistic  for  each  element  in  X 


Function  enoughruns  calls  function  bootstrap_mean  n  times  to  obtain  a  vector  M  of  bootstrapped 
means.  It  then  sorts  the  elements  of  M  in  ascending  order  and  finds  the  (a/2)  percentile  and  (1-a/2) 
percentile  elements  of  M.  If  both  these  %ile  elements  are  within  the  interval  specified  by  function 
interval  above,  (+/- 1 0%  of  the  overall  median)  it  returns  the  value  enoughruns  =  1 ,  otherwise  0. 

enoughruns(n,X,a,y)  :=  j  for  ie  1 .. n 

Xb<-  bootstrap_k(  X ) 

M-*-mean(  Xb) 

med«-median(M) 
perc_intvl«-  perc_interval  (m  ,  a) 
boundvec*-  interval  (med,  y) 
enough  ^<-1  if  boundvec^  percjntvlj 

enough otherwise 
enough ^4- 1  if  perc_intvl2<boundvec2 
enougl^t-O  otherwise 

(stop*- 1)  if  enough= 

stop<-0  otherwise 

perc  JintvU  —  med 

shape* - - - 

med  —  perc_intvlj 

stop 

boundvec ^ 
perc_intvlj 
med 

percjntv^ 
boundvec2 
shape 
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Function  rep_enoughruns  replicates  function  enoughruns  n  times. 


rep_enoughruns  (reps ,  n ,  X ,  a ,  y)  •  - 


for  ie  1 ..  reps 
ER^enoughruns^.X^y)  j 
m<-mean(  ER) 


Function  whatis_enough  determines  how  many  runs  are  enough  using  Candidate  S2 


whatis_enough  (n ,  X ,  a ,  y,  step ,  start) 


r«-rows(X) 


for  i  6  start.,  r 


ER«-0 

if  mod(i,step)sO 

Xi«-  submatrix  (x ,  1 ,  i ,  1 , 2) 
Xiw«-  weightPks(  Xi ) 

ER<-  enoughruns  (n ,  Xi  w ,  a ,  y) 
if  ERj«l 

runs«-i 
mean*— ER^ 

break 

^runs 
jnean 


Function  rep_enoughruns  replicates  function  whatis_enough  n  times. 


rep_whatis_enough  (reps,  n , X , a, y, step, start)  ’ 1 


for  i  6  1 .. 


reps 


ER 


ER 

T 


-  whatis_enough(n ,  X ,  a ,  y,  step,  start) 


Function  pred_std  takes  a  vector  X  and  its  counterpart  Xw  that  is  equally  weighted  according  to 
Firings,  the  significance  level  a  and  the  relative  error  y  and  returns  a  prediction  of  the  number  of 
runs  to  go  using  Candidate  P2,  using  the  bootstrap  estimate  of  standard  error 


pred_boot  (x ,  Xw ,  a ,  yadj ,  se) 


n*-rows(Xw) 
avgfirings«-mean  (x^^ 
0«— mean(  Xw) 


for  j  6  1  ..  500 

n_new«-  n  -f-  floor(  j  •  avgfirings ) 

se_new« — —iV** 

V  n_new 

tvaluei-qt/l  —  —  ,n_new~  1 


hl_new«-tvaluese_new 

. *  hi  new 

ratio_new« - =- — 


Irunstogo«-j 
break  if  ratio__new<yadj 
runstogo 
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Function  bootstrap_se  takes  a  vector  X  and  the  number  of  bootstrap  replications  n  and  returns  the 
bootstrap  estimate  of  standard  error 

bootstrap_se(n,X)  for  ie  l..n 

0bj«-bootmean(X) 

sdb<-stdev(0b) 

boot_se4-sdb 

boot_se 


Function  pred_std  takes  a  vector  X  and  its  counterpart  Xw  that  is  equally  weighted  according  to 
Firings,  the  significance  level  a  and  the  relative  error  y  and  returns  a  prediction  of  the  number  of 
runs  to  go  using  Candidate  PI 

pred_std(x,Xw,a,Y)  !=  n«-rows(Xw) 

avgfirings«-mean  (x<2>) 

0<— mean(Xw) 


hl<-  Clhalflength  (Xw,  a) 


ratio* — ~ 
0 


Y 

-yadj* - 

1-bY 

for  j  €  1  ..  500 


runstogo 

avgfmngs 

0 

n 

ratio 

n_new 

ratio„new 
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Function  stopcrit  takes  a  vector 
X,  the  significance  level  a  and 
the  relative  error  y  and 
determines  when  enough  runs 
have  been  made  using 
Candidate  SI. 


stopcrit (x, a, y) 


q4-rows(X) 
Xs«-csort(x,  l) 


max«-Xs 


‘q.l 

0 


max-)-  1 

7adj« - — 

1+-Y 

for  j  e  20..  q 

Xj4- submatrix  (x,  1 ,  j ,  1 , 2) 
Xj  w«-  weightPks(  Xj ) 
n«-rows(Xjw) 


stdev(Xjw)* 


n—  1 


Vn 


if  n>l 


1020 

se-f-2  otherwise 


6j«-mean(Xjw) 

hl4-qt(l  -  — ,n—  1 

2 

CILoWj«-*0j  —  hi 
CIUpj<—  0j  -(-  hi 


ratiOj< 


hi 


stopf- 1  if  ratiOj<Yadj 
stopt-O  otherwise 
STOP  j  1  if  stop—l 

STOP ^  4-0  otherwise 
STOP24-l  if  j>20 
STOP2+-0  otherwise 


if  STOPS 


Ipredj4-0 
break 

predj  4—  pred_std  (xj ,  Xj  w ,  a ,  yadj)  j 


CILoWj 

ei 

CIUPj 

ratioj 

pred 


otherwise 
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Function  stopboot  determines  when  enough  runs  have  been  made  using  the  bootstrap  estimate 
of  standard  error  and  calls  function  pred_boot  after  each  run  to  predict  the  number  of  runs  to  go 
using  the  bootstrap  estimate  of  standard  error 


stopboot(x,b,a,Y,pilotb)  !-  |q«-rows(X) 


|  for  j  e  pilotb..  q 

Xj<—  submatrix  (x,  1  ,j,  1 ,2) 
Xjwf- weightPks(  Xj ) 
n*-rows(Xjw) 

sej «-  bootstrap_se(  b ,  Xj  w)  if  n>l 
1020 

se-4-2  otherwise 


0j*-mean(Xjw) 
hU-qtfl  -  —  ,n—  iVse. 


CILowj<-0j  —  hi 

CIUPj<-0j-fhl 

ratio.*-— 

J  ft 


stop* —  1  if  ratiOj<yadj 

stop  <-0  otherwise 
STOPjf-1  if  stop=l 

STOPj<-0  otherwise 
STOP2*-l  if  j>20 
STOP2<-0  otherwise 


I  if  STOP= 


lpredj+-0 


|  predj «-  pred_boot  {Xj ,  Xj w ,  a ,  yadj ,  scjj  otherwise 
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ANIMATION/GRAPHING  FUNCTIONS: 


Function  animate_bootstrap_1  replicates  the  bootstrap  resampling  (1  pick  at  a  time) 

of  a  vector  X  in  order  to  show  1  rep  of  the  resampling  process  in  the  animated  graphs  below 


animate  ^bootstrap.  1  ( X ,  p ) 


Xp«-  submatrix  (x,  1  ,p,  1 , 


0 


p«-rows(Xp) 

Xs«-sort(Xp) 


max«-XSp 


W-h*"0 


for  i  €  1 ..  p 
for  j  G  1 ..  p 
Aj  4- bootstraps  ( X ) 
gem*— Aj 


Zjf-z 
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Function  animate_bootstrap_k  replicates  the  bootstrap  resampling  process  (k  picks  at  a  time) 
of  a  vector  X  in  order  to  show  n  reps  of  the  resampling  process  in  the  animated  graphs  below, 
a  is  the  significance  level,  7  is  the  relative  error,  and  p  is  the  number  of  data  points  to  consider. 


animate_bootstrap_k (n,X,a,y,p) 


Xp<- submatrix  (x,  1  ,p,  1 ,  l) 

Xs«-sort(Xp) 

max«-Xs 

P 

zmax-f- 
^max-f- 1*“^ 


yadjj< - I— 

1-hY 


Ilf-0 
12  4- 0 
for  ie  l..n 
for  j  e  1 ..  p 

Aj«- bootstraps  Xp) 

’vr’vr1 

As4-sort(  A) 

0b^f-floor(  mean(  A) ) 

0^4-  floor  (mean  (0b)) 

7adj;<-— L. 

1  +Y 

int«-interval(8,y) 

114-Il-t-ebj 

if  i>l 


£ 

m=  1 


eb„ 


11 


0-1) 

qt (l  -  a,i)  se0^ 


otherwise 


se0^f-O 


J  ratio  j  <-0 
gen.<-A 
YYjf-y 

>,rv‘ 

Z-4-Z 

0- <— mean  (0) 
SE0^4-se0 
RATIO-4- ratio 
F  j4— 7adj 


1 

2 
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|  ©B|<e-0b 
gen 
Z 

YY 
se0 
SE0 
0 
0 

ratio 
RATIO 
7adj 

r 

6b 
0B 


Function  range  determines  the  range  of  data  in  the  subvector  of  size  p  of  the  vector  X  of  size  k  so  the 
range  of  the  animated  graphs  can  be  automatically  set 

range(X,p)  r«-p 

Xp«- submatrix (x,  1  ,p,  1 ,  l) 

Xs«-sort(Xp) 

max*-Xsr 

min<-Xsj 

min  \ 
max/ 
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Function 
animate_stopcrit 
animates  the  process 
of  function  stopcrit. 


anim_stopcrit  (x ,  a ,  y) 


q«-rows(X) 

Xs«-csort(x,  l) 
max«-Xs„  , 

q*1 

zmax+- 
for  j  e  1 ..  q 

Xj«- submatrix  (x,  1 ,  j ,  1 1 2) 
Xqw«—  weightPks(  Xj ) 
n«-rows(Xqw) 


stdev(Xqw) 


n—  1 


J  V"n 

>:<— mean(  Xqw) 


if  n>l 


ratio.  < 


qt|  1  -  ^,n-  1 


0; 


CIlow-<-0j  —  qtf  1  -  —  ,n-  1 

J  J  1  2 


‘sei 


CIhi.4 


7adj:4 


1  -FY 

stop<-  1  if  ratiOj<yadjj 
stop«-0  otherwise 

0j<-0 

SEj*-se 
RATIOjf- ratio 
fadjj^-yadj 


XXj«-X 


<1> 


'(xj.i+1)<_z(xj,i+I) 


se 

SE 

e 

0 

ratio 

RATIO 

7adj 

Tadj 

z 

Z 

CIlow 

Clhi 
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0 

4 

0.4 

5 

to  determine  when  enough  BRAWLER 

o 

4 

0.5 

2 

runs  have  been  made) 

0 

3 

1 

1 

0.4 

5 

0 

1 

0 

0 

0.5 

2 

0.33 

3 

0.5 

2 

0 

2 

0 

1 

0.5 

2 

0 

1 

0 

2 

0 

0 

0 

0 

0 

2 

0.25 

4 

0.33 

3 

0.33 

3 

0.5 

2 

1 

1 

0 

0 

0.5 

2 

0.67 

3 

GLOBAL  PARAMETERS: 

0 

2 

0 

0 

0.5 

4 

1 

2 

B  :=501 

0 

3 

0 

0 

0.14 

7 

p  =99 

1 

2 

0 

0 

0 

4 

1 

1 

q  =122 

0.5 

4 

0.5 

2 

0 

1 

0.5 

4 

0.33 

3 

a  =.10 

1 

2 

1 

1 

0.67 

3 

0 

3 

0.5 

2 

0 

7 

y  :-.10 

1 

2 

1 

2 

0 

5 

0 

1 

0.14 

7 

0 

0 

0 

3 

0.17 

6 

0.33 

3 

0 

1 

0.33 

3 

0.67 

3 

0 

2 

0.2 

5 

0 

1 

0.5 

4 

0 

0 

— 

- 
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Pki  :  = 


0 

0 

1 

0 

0 

0.25 

0 

0.5 

1 

0 

0 

0 

0.25 

1 

1 

0 

1 

1 

0.5 

0 

0 

0 

1 

0.33 

1 

1 

1 

0.17 

0.33 

0.14 

1 

1 

0 

0 

0 

0.5 

0.25 

0 

0 

0 

0 

1 

0 


1 

0 

2 

0 

2 

4 

4 

2 

1 

2 

0 

1 

4 

3 

1 

1 

1 

1 

2 

2 

1 

2 

1 

3 

2 

1 

1 

6 

3 
7 
1 
1 
0 
0 

4 
6 

4 

5 
1 
2 
2 
1 
0 


Y  is  the  output  data  Adjusted 

Kills  per  Firing  and  Firings  from  BRAWLER 


Datal  :  -( stack(  Pk  1 ,  Pk2 )  )•  1 00 
Data2  1  =stack(  FI ,  F2 ) 

Y  !  =augment(  Datal ,  Data2) 

Y  is  the  BRAWLER  data  matrix 
of  AKPPs  and  Firings. 

Y  must  be  made  of  integers 
only,  therefore  the  AKPF's 
are  multiplied  by  100 
because  Mathcad 
animated  simulations  only 
work  with  integer  data 


till 

SI  40 

5 

m  33 

3 

SI  50 

4 

50 

4 

fH  20 

5 

Is  25 

4 

ji!  100 

1 

||  40 

5 

ifeso 

2 

§0: 100 

1 

11  0 

1 

12  50 

2 

&  50 

2 

i§o 

1 

lo 

1 

k  :  -rows(  Y) 
k  =  200 

k  is  the  #  of 
data  points  (runs) 


u 

l 

0.5 

2 

0 

1 

0.75 

4 

0.17 

6 

0.33 

3 

0.5 

2 

0.25 

4 

0 

1 

0 

F2  := 

0 

0.33 

3 

0.33 

3 

0 

0 

0.25 

4 

0.33 

3 

0 

2 

0 

3 

0.5 

2 

1 

2 

0.67 

3 

0.4 

5 

0.25 

8 

0 

1 

0 

4 

0.5 

2 

0.75 

4 

0.33 

6 

0.5 

4 

0 

o 

1 

1 

0 

2 

0.2 

5 

0.33 

3 

0 

1 

0 

1 

0.33 

3 

0.5 

2 

0.43 

7 

0.2 

5 

1 

2 

0 

3 

0.33 

3 

0 

2 

0 
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Yp  ! -submatrix (y,  1  ,p,  1 ,2) 


Yq  : -submatrix  (y,  1  ,q,  1 ,2) 


Ys  :-csort(Yp,  l) 


Yw  :-weightPks(  Y) 
sumfirings(  Y)  =  506 

mean(Yw)  =32.178 


Ypw  ’=weightPks(  Yp) 
sumfirings(Yp)  =229 
mean(Ypw)  =34.048 


Yqw  :-weightPks(  Yq) 
sumfirings(  Yq)  =  287 

mean(Yqw)  =32.387 


Yqw 
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Clvecl  :=boot_CIs(Y,B,a,l) 


meanvecl  ! -Clvecl 


Cllowvecl  :-CIvecl0 


meanvecl  = 


1  1  1 

|®40 

fc  37.375 

■  41.583 

M  43.688 

j|  38.19 

It  35.8 

It  37.923 

Hi  38.484 

Si  39.212 

B  41.088 

I®  39.771 

If 40.622 

P  40.897 

P  40.05 

K  39.171 

enoughruns(B ,  Ypw,a,y)  = 


Clupvecl  :-CIvecl3 


v  38.5 


4  40.75 


Cllowvecl ! 


0 

30.628 

30.843 

34.031 

37.817 

37.434 

1.188 


II  40 


11  35.625 


IS  33.762 


■  32.16 


Clupvecl  = 


P I  33.269 


13  34.484 


!  35.061 


|  36.5 


m  35.057 


l  35.459 


3  36.769 


i  35.25 


■  34.024 


- 

p40 

|§j  39.125 

ll  44.917 

U  46.438 

fe'-:  42.19 

I!  40.04 

■  43.808 

43.258 

IP  44.03 

10  46.853 

ft  45.257 

■  45.432 

B  45.718 

■  45.65 

i5: 44.171 

boot_Z(Ypw,B)  = 


fe  j 

in 

2^| 

lOit 

un 

5; 

if! 

P 

28 

.62 

1.941 

-2.797 

34 

048 

2.072 

39 

.843 

ip! 

28 

.803 

1.939 

-2.705 

34 

048 

2.072 

39 

653 

Bj 

29 

.013 

1.981 

-2.542 

34 

048 

2.072 

39 

314 

m 

28 

983 

2 

-2.533 

34 

048 

2.072 

39 

297 

S 

29 

218 

1.938 

-2.492 

34 

048 

2.072 

39 

212 

29 

437 

1.903 

-2.424 

34 

048 

2.072 

39 

07 

rj 

29 

38 

1.988 

-2.348 

34 

048 

2.072 

38 

913 

29. 

.738 

1.841 

-2.341 

34. 

.048 

2.072 

38. 

.899 

29, 

.886 

1.806 

-2.304 

34. 

.048 

2.072 

38, 

.823 

!i 

30, 

.109 

1.751 

-2.25 

34, 

.048 

2.072 

38 

.71 

■ 

29 

.694 

1.941 

-2.244 

34, 

.048 

2.072 

38 

.697 

m 

29 

.703 

1.981 

-2.193 

34, 

.048 

2.072 

38, 

.593 

w. 

30 

1.848 

-2.19 

34, 

.048 

2.072 

38 

.586 

n 

29 

.738 

1.976 

-2.182 

34. 

.048 

2.072 

38 

.569 

f5" 

29 

.913 

1.929 

-2.144 

34, 

.048 

2.072 

38 

.49 
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rep_enoughruns(l ,  B ,  Ypw ,  a ,  y)  =  1 


whatis_enough(B ,  Y,  a,y,  1 , 1 80) 


180  \ 
32.639  / 


rep_whatis_enough(l  ,B,  Y.a.y,  1 , 180)  =(  180  32.704  ) 


stopcrit(Y,a,y) 


122 

29.448184 

32.38676 

35.325335 

0.090734 

{-7681 70609, 1 073080493 } 


res  :-animate_bootstrap_l(Yw,p) 

r  1  =Ys  j ,  i  + 1  ••  Ysp,  1 4- 1 


gen^ -res  j  Z^-reSj 
i:=res3 


Animated  Bootstrap  Resampling  Animated  Bootstrap  Resampling  Histogram 

choose  k  out  of  k  with  replacement  (Cumulative)  choose  k  out  of  k  with 

(1  frame  =  1  choose)  replacement  (1  frame  =  1  choose) 


genFRAME+1,r  =  values  of  bootstrapped  data  points 


r  =  value  of  bootstrapped  data  point 
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j  ”1  ..p 

res  '  -animate_bootstrap_k  (b  ,  Y,  a ,  y,  p) 

b:  =  l..B 

gen  ! -res j  Z:-res2 

YY  '  -res-j  se0  !  -res^ 

SE0  :  -res^ 

0  -res^ 

Ti  ~resjQ  r:=resn 

0b:=res^2  0B  -res^ 

0  !  -res^ 

ratio  '  “resg 

range(X,p)  :  = 

r«-P 

Xp<- submatrix (x,  1  ,p,  1 ,  l) 

Xs«-sort(Xp) 

max4-Xsr 

min<-Xsj 

/  o  \ 

range(  Y,p)  = 

\  100/ 

RATIO !  -resg 

/min 

\max 


r : -range(  Y,  p )  j  +  1 ..  range(  Y,  p)2  -f- 1 


Animated  Resampling 

B  bootstrap  samples  (1  frame  =  1  bootstrap 
sample  for  each  of  B  bootstrap  reps) 


T 

• 

•  • 

• 

• 

•t •#  mm 

•  • 

•  - 

• 

m  •  • 

»  •  • 

•  •  • 

• 

-J _ 

•  •  • 

• 

•  • 

• 

•  • 

m 

pi  ■■■  ■  » — — ..  J 

0  50  100 


Cumulative  Animated 
Bootstrap  Sample  Histogram 

B  bootstrap  samples  (1  frame  =  1  bootstrap 
sample  for  each  of  B  bootstrap  reps) 


Animated  Bootstrap  Means  (Cumulative) 

B  bootstrap  means  (1  frame  =  1  median 
for  each  of  B  bootstrap  reps) 
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Cllowvec  •  -anim_stopcrit(Y, a, y) j  j 


mean vec  ‘  -  anim_stopcrit  (y  ,  a ,  y)  3 


Clupvec  :-anim_stopcrit(Y,a,  y)  12 


meanvec 


Si  40 

|§j  40 

1  O 

R  37.375 

nil 

H  34.781 

|f  39.969 

H  41.583 

H  37.871 

K  45.295 

1 43.688 

S  40.531 

H  46.844 

111  38.048 

@!  33.43 

B  42.665 

1^;:; 

35.96  Cllowvec  = 

Jit  31.78  Clupvec  = 

1 40.14 

Hi  38.423 

H  32.552 

1 44.294 

H  38.677 

u> 

w 

bo 

!§  43.542 

1  39.364 

Is  34.739 

fe  43.989 

41.147 

H  35.719 

fe  46.575 

H  39.971 

|1  34.334 

||  45.608 

H?  40.514 

H  35.159 

If  45.868 

■  41 

11  35.901 

ji  46.099 

B  39.975 

H  34.711 

P  45.239 

Hi  39 

ft  33.608 

IS  44.392 

j  :  1  •*  q 


res  •  ~anim_stopcrit  (y  ,  a ,  y) 


=Ysl,l”Ysp,l 


b:=l..B 


se  ^-resj 
yadj  1  -res  ^ 


SE:-res2  0  ^-res^  0  ’-res4  ratio  :-res5 
T adj  1  “resg  z  ■  -res^  Z  i  =res  jq 


RATIO  !  -res^ 


r  : -range(  Y,q)  j -j- 1 ..  range(  Y,  q)2-f- 1 


Animated  Data  from  Runs 


Animated  Data  from  Runs  Histogram 


JL 

100 


150 
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Appendix  B.  enoughrunsl 7 


program  enoughruns 

C***  THIS  ROUTINE  IS  WRITTEN  AT  RAND 

C# ABSTRACT  Program  to  test  whether  enough  tbmain  runs  have  been  made 
C#PURPOSE  Find  the  mean  and  confidence  interval  for  exchange  ratio 
C# AUTHOR  Carol  Johnson  (bootstrap  method  from  Tom  Lucas,  RAND) 

C#TYPE  Postprocessing 

C# PARAMETER  DESCRIPTIONS: 

CIN  inarg(l)  CHAR  -  (optional)  ' eratio'  or  'red'  or  'blue' 

C# TECHNICAL  DESCRIPTION: 

C  Program  enoughruns  tells  whether  "enough"  tbmain  runs  have  been  made 
C  Used  by  the  fly  script. 

C  Input  file:  fly. TOTALS,  with  lines  like: 

C : 0  1  BLUE  AC  +  1  RED  AC  =  2  TOTAL.  0  BLUE  AC  ALIVE  +  . 

C  Output  file:  fly. PERMEAS 

C  stdout :  "yes,  _  runs  are  enough"  or  "no,  _  runs  are  not  enough" 

C#VARIABLE  DESCRIPTIONS: 

C  MAXJRUNS  PAR  -  maximum  number  of  lines  in  fly. TOTALS  to  process 
C  iargc  INT  FUNCTION  -  number  of  arguments  on  the  command  line 
C  nargs  INT  -  number  of  arguments  on  the  command  line 

C  indexpm  INT  -  1,  2,  or  3  for  3  argument  possibilities 

C  num_runs  INT  -  number  of  lines  in  fly. TOTALS 

C  bac  INT  -  number  of  blue  a/c 

C  balive  INT  -  number  of  blue  a/c  alive  at  the  end  of  this  run 

C  bkrun  INT  array  -  number  of  blue  a/c  killed,  in  each  run 

C  rac  INT  -  number  of  red  a/c 

C  ralive  INT  -  number  of  red  a/c  alive  at  the  end  of  this  run 

C  rkrun  INT  array  -  number  of  red  a/c  killed,  in  each  run 

C  pmmean  REAL  -  performance  measure's  mean 

C  pramin  REAL  -  performance  measure's  least  value 

C  pmlow  REAL  -  performance  measure's  5th  percentile 

C  pmmedian  REAL  -  performance  measure's  median 

C  pmhigh  REAL  -  performance  measure's  95th  percentile 

C  pmmax  REAL  -  performance  measure's  greatest  value 

C  inarg(l)  CHAR  -  argument,  or  default  if  no  argument 
C  linel  CHAR  -  a  line  of  the  fly. TOTALS  file,  as  read  in 
C  answer  CHAR  -  'yes'  or  ' no ’ ,  enough  runs 
C#### 

C#AUDIT 

C  BY  tacb  IN  v6.15  tacb. cruise  ON  95/06/16  18:09:27  to.n.rou 
C  BY  tacb  IN  v6.15  tacb. cruise  ON  95/06/16  17:57:04  MODIFIED 
C  Putting  illegal  unit  -1  into  a  "variable",  to  satisfy  H-P 
C  BY  tacb  IN  v6.15  skip.cnj  ON  94/12/28  13:53:03  to.n.rou 

C  BY  tacb  IN  v6.15  skip.cnj  ON  94/12/21  13:50:59  MODIFIED 

C  Split  out  subroutine  bootstrap,  to  call  from  analyze  pgm. 

C  BY  tacb  IN  v6.15  skip.cnj  ON  94/12/21  13:47:55  to.n.rou 

C  BY  tacb  IN  v6 . 15  skip.cnj  ON  94/12/20  13:43:56  MODIFIED 

C  Add  argument  to  choose  measure  of  performance 
C  BY  tacb  IN  v6.15  tacb. rand  ON  94/10/17  14:50:31  to.n.rou 

C  BY  tacb  IN  v6.15  tacb. rand  ON  94/10/10  16:45:37  MODIFIED 

C  Change  to  bootstrap  method  a  la  Tom  Lucas,  RAND. 

C  Loop  to  find  minimum  number  of  "enough"  runs. 

C  Add  field  "#Runs"  to  EXCHANGE  RATIO  STATISTICS 
C  BY  tacb  IN  v6 . 15  tacb. rand  ON  94/09/20  17:01:56  to.n.rou 

C  BY  tacb  IN  v6.15  tacb. rand  ON  94/09/20  14:21:39  MODIFIED 

C  First  test  was  of  runs  1-11;  now  is  of  runs  10-20 
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C  BY  tacb  IN  v6.15  tacb. analyze  ON  94/09/13  11:56:46  to.n.rou 
C  BY  tacb  IN  v6.15  tacb. analyze  ON  94/09/13  11:49:26  MODIFIED 
C  Expect  every  line  to  be  a  TOTAL,  line 

C  BY  tacb  IN  v6.15  tacb. analyze  ON  94/09/13  11:18:54  to.n.rou 
C#### 

integer  MAX_RUNS 
parameter (MAX_RUNS=1001) 
integer  iargc ,  nargs , indexpm 
integer  num_runs , i , minus 1 
integer  bac ,  balive ,  bkrun  (MAXJOJNS) 
integer  rac , r alive , rkrun (MAX_RUNS) 
real  pmmean , pmmin , pmlow, pmmedian , pmhigh , pmmax 
character *70  inarg(l) 
character*80  linel 
character* 3  answer 
*  ENDDEC 

--GET  ARGUMENTS 
minusl  =  -1 
nargs=iargc ( ) 
if  (nargs.eq.O)  then 
inarg(l)  =  'eratio' 
indexpm  =  1 
else 

call  getargd,  inarg(l)  ) 
if  (inarg(l) .eq. 'eratio’ )  then 
indexpm  =  1 

elseif  (inarg(l) .eq. 'red' )  then 
indexpm  =  2 

elseif  (inarg(l) .eq. 'blue ' )  then 
indexpm  =  3 
else 

write (6,*)  'NOTE  analyze  gets  unknown  argument  '//inarg(l) 
write (6,*)  '  so  writing  illegal  unit  -1  to  abort' 

write (minusl, *)  '  so  writing  illegal  unit  -1  to  abort' 

endif 
endif 

C  --INITIALIZE  FILES  AND  COUNTS  AND  SUCH 

open (unit=5 , f ile= ' fly . TOTALS ' , status- ' old' , err=2000 ) 
open (unit=ll , f ile= ' fly. PERMEAS ' , status= ' unknown ' , err=3000 ) 

C  — READ  NUMBERS  OF  BLUE  AND  RED  AIRCRAFT  KILLED  IN  SEVERAL  RUNS 

do  100  num_runs=l , MAX_RUNS 
read ( 5 , fmt= ' (a) ' , end=200 , err=2000)  linel 
if  (linel (40:45)  .eq.  'TOTAL.'  )  then 
C  --  Find  number  of  aircraft  killed  in  this  run 

read (linel, ' (  tl2 , i2 , t25 , i2 , t49 , i2 , 17x, i2 ) ' ) 

>  bac,  rac,  balive,  ralive 
bkrun (num_runs)  =  bac  -  balive 
rkrun (num_runs)  =  rac  -  ralive 

else 

goto  2100 
endif 

100  continue 
200  continue 

num__r un  s = num_r un  s  - 1 

C  --BOOTSTRAP,  IF  A  MINIMUM  NUMBER  OF  RUNS 

if  (num_runs . ge . 20)  then 

call  bootstrap ( ’ enoughruns ' , indexpm, 11 , num_runs , bkrun, rkrun, 20 , 

>  pmmean ,  pmmin ,  pmlow ,  pmmedian ,  pmhigh ,  pmmax ,  answer ) 
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/  i3 , 


runs  are  not  enough' ' ) ' ) 


else 

answer  =  ' no ' 
write (11,*)  '  ' 

write ( 11 , fmt= ' ( ' 'no, 

>  num_runs 

write ( 6 , fmt= ' ( 1 ' no ,  ' ' , i3 , ' '  runs  are  not  enough ' ' ) ' ) 

>  num_runs 
endif 

stop 

2000  write  (6,2001) 

2001  format  ( ' enoughruns :  File  fly. TOTALS  is  not  readable') 

goto  9000 

2100  write  (6,2101)  linel 

2101  format  ('enoughruns:  Expected  TOTAL  line  but  found' /a) 

goto  9000 

3000  write  (6,3001) 

3001  format  ('enoughruns:  File  fly. PERMEAS  is  not  writable') 

9000  write (6,  *)  '  so  writing  illegal  unit  -1  to  abort' 

write (minusl, *)  '  so  writing  illegal  unit  -1  to  abort' 

end 
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subroutine  bootstrap (caller , ipm, iunit , num_runs , bkrun, rkrun, 

>  min_run ,  pmmean ,  pmmin ,  pmlow,  pmmedian ,  pmhigh ,  pmmax ,  answer ) 

C***  THIS  SUBROUTINE  IS  WRITTEN  AT  RAND 

C# ABSTRACT  Get  statistics  on  exchange  ratio  or  other  performance  measure 
C# PURPOSE  Help  determine  whether  enough  runs  have  been  made 
C# AUTHOR  Carol  N.  Johnson  {bootstrap  method  from  Tom  Lucas,  RAND) 

C#TYPE  Post-processing 

C # PARAMETER  DESCRIPTIONS: 

CIN  caller  CHAR  -  name  of  the  calling  program  (analyze  or  enoughruns) 
CIN  ipm  INT  -  performance  measure  option: 

C  1  =  eratio; 

C  2  =  red  a/c  killed; 

C  3  =  blue  a/c  killed 

CIN  iunit  INT  -  Fortran  logical  unit  at  which  to  write 

CIN  num_runs  INT  -  maximum  number  of  runs  to  bootstrap 

CIN  bkrun (num„runs)  INT  -  number  of  blue  a/c  killed,  in  each  run 

CIN  rkrun(num_runs)  INT  -  number  of  red  a/c  killed,  in  each  run 

CIN  min__run  INT  -  minimum  number  of  runs  to  bootstrap 

COUT  pmmean  REAL  -  performance  measure’s  mean 

COUT  pmmin  REAL  -  performance  measure's  least  value 

COUT  pmlow  REAL  -  performance  measure's  5th  percentile 

COUT  pmmedian  REAL  -  performance  measure's  median 

COUT  pmhigh  REAL  -  performance  measure's  95th  percentile 

COUT  pmmax  REAL  -  performance  measure's  greatest  value 

COUT  answer  CHAR  -  ' yes '  or  ' no 1 ,  enough  runs 

C#TECHNICAL  DESCRIPTION: 

C  This  is  a  poor  man's  way  of  making  N  runs  501  times: 

C  Given  numbers  of  red  and  blue  aircraft  killed  in  each  of  N  runs: 

C  do  501  (or  2001)  times: 

C  randomly  draw  N  of  N  runs 

C  calculate  exchange  ratio  (or  use  other  measure  of  performance) 

C  sort  501  exchange  ratios 
C  eratio (251)  is  the  median 

C  90%  of  the  501  eratios  are  between  eratio (26)  and  eratio (476) 

C  when  90%  of  eratios  are  within  the  median  plus  or  minus  10%, 

C  N  runs  are  enough 
C# VARIABLE  DESCRIPTIONS: 

C  ITERATIONS  PAR  -  number  of  times  to  simulate  making  num_runs  runs 
C  iseq  INT  -  random  number  seed 

C  i  INT  -  index 

C  now_runs  INT  -  index  for  runs 

C  ib  INT  -  index  for  bootstrap  iterations 

C  index  INT  -  randomly- chosen  run 

C  bksum  INT  -  sum  of  blue  a/c  killed  in  num_runs  random  runs 

C  rksum  INT  -  sum  of  red  a/c  killed  in  num_runs  random  runs 

C  pmsum  REAL  -  sum  of  performance  measure  in  num_runs  runs 

C  permeas  REAL  -  performance  measure  in  each  bootstrap  iteration 

C  anti  LOGICAL  -  .true,  or  .false.,  antithetical  random  numbers 

C#### 

C# AUDIT 

C  BY  tacb  IN  v6.15  tacb. analyze  ON  95/10/25  18:42:21  to.n.rou 

C  BY  tacb  IN  v6.15  tacb. analyze  ON  95/10/25  17:25:41  MODIFIED 

C  Avoid  NaN  when  median  is  zero. 

C  BY  tacb  IN  v6.15  tacb. analyze  ON  95/09/18  10:19:11  to.n.rou 

C  BY  tacb  IN  v6.15  tacb. analyze  ON  95/09/15  14:06:57  MODIFIED 

C  For  analyze  program  and  ipm=6,  print  combined  report  for  ipm=l+2+3 

C  BY  tacb  IN  v6.15  tacb. analyze  ON  95/08/02  13:47:25  to.n.rou 

C  BY  tacb  IN  v6.15  tacb. analyze  ON  95/08/02  11:00:27  MODIFIED 

C  Print  which  kind  of  performance  measure  is  being  used,  in  analyze  pgm 

C  Change  indexpm  to  ipm 
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C  BY  tacb  IN  v6.15  skip.cnj  ON  95/01/31  16:40:32  to.n.rou 
C  BY  tacb  IN  v6 . 15  skip.cnj  ON  95/01/31  16:31:14  MODIFIED 
C  If  caller  is  analyze,  say  no,  not  enough  if  fewer  than  20  runs 
C  BY  tacb  IN  v6.15  skip.cnj  ON  94/12/28  14:09:50  to.n.rou 
C#### 

C  Args : 

character* ( * )  caller 
integer  ipm, iunit 

integer  num_runs,  bkrun ( num_runs ) ,  rkrun ( num_runs ) ,  min_run 
real  pirimean,  pmmin,  pmlow,  pmmedian,  pmhigh,  pmmax 
character* ( * )  answer 
C  Functions : 

real  rnsq 
C  Locals : 

integer  ITERATIONS 
parameter (ITERATIONS=  501) 

CUT  parameter ( ITERATIONS=2001) 
integer  iseq 

integer  i , now_runs , ib , index 
integer  bksum  , rksum 
real  pmsum 

real  permeas ( ITERATIONS) 
logical  anti 
C*ENDDEC 

iseq  =  0 
anti  =  .false, 
call  ranini (iseq, anti) 
do  7  00  now_runs=num__runs  , min_run,  -1 
C  --ITERATE,  CALCULATING  PERFORMANCE  MEASURE 

pmsum  =  0 . 

do  500  ib=l, ITERATIONS 
bksum=0 
rksum- 0 

do  220  i  =  1,  now_runs 

index  =  NINT (rnsq ( iseq) *real (now_runs )+. 500001 ) 
bksum=bksum  +bkrun( index) 
rksum=rksum  +rkrun { index) 

220  continue 

goto  (310,320,330),  ipm 
310  continue 

if  (bksum. ne . 0 . )  then 

permeas (ib)  =  float (rksum) / float (bksum) 
else 

permeas (ib)  =  l.E+20 
endif 
goto  390 

320  continue 

permeas (ib)  =  float (rksum) /float (now_runs) 
goto  390 

330  continue 

permeas  (ib)  =  float  (bksum)  /  float  (now__runs) 
goto  390 

390  continue 

pmsum  =  pmsum  +permeas ( ib) 

500  continue 

pmmean  =  pmsum/ real  ( ITERATIONS) 

C  --SORT  PERFORMANCE  MEASURES  ASCENDING,  IN  PLACE 

call  sort (permeas, ITERATIONS) 
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pmmin  =  permeas ( 1 ) 

pmmedian  =  permeas (( ITERATIONS+1 ) /2 ) 
pmlow  =  permeas (ITERATIONS/2 0) 
pmhigh=  permeas ( ITERATIONS  -ITERATIONS/2 0 ) 
pmmax  =  permeas ( ITERATIONS) 

C  --DECIDE  WHETHER  ENOUGH  RUNS 

if  (caller  .eq.  'analyze' 

>  . and .  now_r uns  .It.  20)  then 

answer  =  1  no ' 

elseif  (pmlow  .ge.  pmmedian*. 90 

>  .and.pmhigh  .le.  pmmedian* 1 . 10 )  then 
answer  =  'yes' 

else 

answer  =  'no 1 
endif 

if  (caller  .eq.  'analyze')  then 

call  writesix  ( ipm, iunit , now_runs , 

>  pmmean, pmlow, pmmedian, pmhigh, permeas , ITERATIONS , answer) 
else 

call  writenow  ( ipm, iunit , now_runs , 

>  pmmean, pmlow, pmmedian, pmhigh, permeas , ITERATIONS , answer) 
endif 

if  (answer  .eq.  'no')  goto  800 
700  continue 

800  continue 

return 
end 
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C - 

subroutine  writenow  ( ipm, iunit , now_runs , 

>  pmmean,  pmlow,  pmmedian,  pmhigh,  permeas,  ITERATIONS,  answer) 

C# PARAMETERS 

integer  ipm, iunit , now_runs , ITERATIONS 

real  pmmean , pmlow , pmmedian , pmhigh , permeas ( ITERATIONS ) 
character* ( * )  answer 
C# LOCAL 

real  FRACTINC , FRACTMAX1 

C  --  FRACTINC,  FRACTMAX1  depend  on  ITERATIONS 

parameter (FRACTINC= . 01 ,  FRACTMAX1= . 905 ) 

CUT  parameter (FRACTINC= . 005 , FRACTMAX1= . 9025 ) 
integer  i , ib , interval s , tvalpop , sumpop , j 
real  fractmin, fractmax, tvalmax 
character* 14  pmname(3) 
data  pmname/ 'EXCHANGE  RATIO', 

>  ' REDS  KILLED ' , 

>  'BLUES  KILLED'/ 

C#ENDDEC 

C  --WRITE  CONFIDENCE  INTERVAL  AND  MEDIAN  OF  PERFORMANCE  MEASURES 


4010  formate  BOOTSTRAP  STATISTICS  —  ',a// 

>  1  < - Extrema  of  90% - >  ’  / 

>  '  5th  95th  '/ 

>  '  #Runs  Mean  Min  %ile  Median  %ile  Max’/ 

>  .  -  -  -  -  -  -  - «/ 

>  i6,6f8.2,’  Value') 

4011  format (14x, 5f 8 . 3 , '  Fraction  of  Median’) 


write ( iunit , 4010 )  pmname ( ipm) , now_runs , pmmean, 

>  permeas (1) , pmlow, pmmedian, pmhigh, permeas ( ITERATIONS) 
if  (pmmedian  .ne.  0)  then 

write (iunit , 4011) 

>  permeas  (1)  /pmmedian,  pmlow /pmmedian,  pmmedian /pmmedian, 

>  pmhigh /pmmedian , permeas ( ITERATIONS ) /pmmedian 
endif 

C  --WRITE  HISTOGRAM  OF  PERFORMANCE  MEASURES 

4020  format (/t9, ’BOOTSTRAP  INTERVALS  —  ’,a// 

>  1  Fraction  of  Median’/ 

>  ’  Min.  <  Middle  <  Max.  Count’) 

4022  format ( f 7 . 3 , '  < ' , f 7 . 2 , '  < ’ , f 8 . 3 , 2x, 200al/ ( 28x, 200al ) ) 

if  (pmmedian  .ne.  0)  then 
write (iunit , 4020)  pmname (ipm) 
fractmin  =  permeas ( 1 ) /pmmedian 

intervals  =  2 *( 1 . 00-FRACTMAX1) /FRACTINC  +2.0005 

fractmax  =  FRACTMAX 1 

tvalmax  =  pmmedian  * fractmax 

ib  =  1 

sumpop  =  0 

do  600  i=l , intervals 
tvalpop  =  0 

580  if  (ib  .le.  ITERATIONS)  then 

if  (permeas (ib)  .le.  tvalmax)  then 
tvalpop  =  tvalpop  +1 
ib  =  ib  +1 
goto  580 
endif 
endi  f 

write (iunit , 4022 )  fractmin, ( fractmin+ fractmax) *. 5 , fractmax, 

>  ( ' + ' , j=l, tvalpop) 
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sumpop  =  sumpop  +tvalpop 
fractmin  =  fractmax 
if  (i.lt. ( intervals-1 ) )  then 

fractmax  =  fractmax  +FRACTINC 
tvalmax  =  pmmedian  * fractmax 
else 

fractmax  =  permeas ( ITERATIONS) /pmmedian 
tvalmax  =  permeas (ITERATIONS) 
end  if 

600  continue 

if  (sumpop  .ne.  ITERATIONS)  then 

write (6, ' ( ' ' WRITENOW  error:  sumpop  .ne.  ITERATIONS' ' ) ' ) 
stop  1 
endif 
endif 

C  --PRINT  WHETHER  ENOUGH  RUNS 

if  (answer .  eq. 'yes ' )  then 

write ( iunit , fmt= ' ( / ' ' yes ,  '  ' , i3 ,  '  '  runs  are  enough ''/)') 

>  now_runs 

write ( 6 , fmt= '  ( '  ' yes  ,  '  ' , i3 ,  '  '  runs  are  enough '  ' )  ' ) 

>  now_runs 
else 

write ( iunit , fmt= 1 (/' 'no,  '  ' , i3 ,  '  '  runs  are  not  enough ' ' ) ' ) 

>  now_runs 

write ( 6 , fmt= ' ( ’ 'no,  ' ' , i3 , ' '  runs  are  not  enough ' ’ ) ' ) 

>  now_runs 
endif 
return 

end 
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C - 

subroutine  writesix  ( ipm, iunit , now^runs , 

>  pminean,  pmlow,  pmmedian, pmhigh, permeas ,  ITERATIONS ,  answer) 
C# PARAMETERS 

integer  ipm, iunit , now_runs , ITERATIONS 

real  pmmean, pmlow, pmmedian, pmhigh, permeas (ITERATIONS) 
character* ( * )  answer 
C# LOCAL 

character* 14  pmname(3) 

data  pmname /' Exchange  Ratio', 


>  'Reds  Killed', 

>  ’Blues  Killed’/ 

C#ENDDEC 

C  --WRITE  CONFIDENCE  INTERVALS  BY  BOOTSTRAPPING 

4010  format  ( 

>  *  CONFIDENCE  INTERVALS  BY  BOOTSTRAPPING'// 

>  '  < - Extrema  of  90% - >  '/ 

>  '  Enough  5th  95th  '/ 

>  '  Runs?  Mean  Min  %ile  Median  %ile  Max') 

4012  format ( 

>  ■  -  -  -  -  -  -  - 

>  '  - ') 

4014  format (2x,a3,2x,f8.2,5f8.2,2x,a) 

4015  format(7x,  8x,  5f8.3,'  Fraction  of  Median') 


i f  ( ipm  . eq .  1 )  then 
write (iunit, *)  *  ' 

write (iunit, 4010) 
write ( iunit , 4012 ) 
endif 

write ( iunit , 4014 )  answer , pmme an, 

>  permeas(l) , pmlow, pmmedian, pmhigh, permeas ( ITERATIONS) , 

>  pmname (ipm) 

if  (pmmedian  . ne.  0)  then 
write ( iunit , 4015 ) 

>  permeas ( 1 ) /pmmedian , pmlow/ pmmedian , pmmedian /pmmedian , 

>  pmhigh /pmmedian , permeas ( ITERATIONS ) /pmmedian 
endif 

write ( iunit, 4012 ) 

return 

end 

c - 
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#  Purpose:  STATS  script  to  post-process  Brawler  runs. 

#  Author (s):  LAWRENCE  J.  TARANTO 

#  JIM  WILLIAMS  (SVERDRUP) 

#  ASC/XRA,  BUILDING  16 

#  2275  D  STREET  SUITE  10 

#  WRIGHT -PATTERSON  AFB  OH  45433-7227 

#  Date:  14  OCTOBER  1997 

#  Notes:  Compatible  with  Brawler  V6.3 

# 

#  Determine  terminal  connection 

set  Terminal=' tty  |  sed  'sr'Vdev/::'  ' 

if  (  "X$Terminal "  ==  "X"  ||  "X$Terminal"  ==  "Xnot  a  tty")  set  Terminal=null 

if  ("$Terminal"  ==  "null")  then 
set  Terminal=" /dev/null " 
else 

set  Terminal= " /dev/ $Terminal " 
endif 

# 

#  Set  local  variables 
set  No=0 

set  RunPath="$ STATS" 
set  Yes=l 

# 

#  Initialize 

set  OverWrite=$Yes 
set  PreProcess=$Yes 
set  Prompts=$No 

# 

#  Process  arguments 
if  ("$1"  ==  " " )  then 

set  ManyDirList= " " 
else  if  ("$1"  ==  "batch")  then 
set  ManyDirList= " " 
set  Prompts=$No 
else 

set  ManyDirList=$l 
endif 

# 

#  Set  date  and  time 

set  CurrentDate=" 'date  1 +%d%b%y'x" 
set  CurrentTime= " ' date  ‘ +%H%M . %S 1 ' " 

# 

echo  "  "  >  $Terminal 

echo  "Brawler  Post-Processor"  >  $Terminal 

echo  "  Execution  Date/Time:  $CurrentDate  /  $CurrentTime"  >  $Terminal 
echo  "  "  >  $Terminal 
# 

#  Determine  if  any  "many"  directories  exist 
if  ( " $ManyDirList"  ==  "")  then 
Is  -d  many*  >&  /dev/null 
if  ($status  ! =  0)  then 

#  Determine  if  any  LOG  and/or  database  files  exist  in  current  directory 

Is  LOG.*  >&  /dev/null 

set  LOGstatus=$status 

Is  database  >&  /dev/null 

set  DBstatus=$status 

if  ( $LOGstatus  ==  0  |  $DBstatus  ==  0)  then 
set  ManyDirList= " . " 
else 

echo  "  " 
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echo  ’ERROR:  No  "many"  directories,  "LOG”  files,  or  "database”  files 

found. . . 1 

echo  "  " 
exit 
endif 
else 

set  ManyDirList=' Is  -d  many*' 
endif 
endif 
# 

if  ($ Prompts  ==  $Yes)  then 

csh  -f  $RunPath/yesno . csh  "If  Database  Files  Exists,  Overwrite"  $Yes  3 
$Terminal 

set  OverWrite=$status 
# 

if  ($OverWrite  ==  $Yes)  then 

csh  -f  $RunPath/yesno . csh  "Pre-Process  LOG  Files”  $Yes  3  $Terminal 
set  PreProcess=$status 
else 

set  PreProcess=$No 
endif 
endif 
# 

set  CurrentDir='pwd' 
foreach  ManyDir  ($ManyDirList) 
cd  $CurrentDir/$ManyDir 
echo  "  *  >  $Terminal 

echo  ‘  Processing  In  Directory  " ’ 'basename  $cwd' ’ . . . 1  >  $Terminal 
# 

Is  LOG.*  >&  /dev/null 
if  ( $status  ! =  0)  then 
echo  "  " 

echo  "ERROR:  No  LOG  files  found. . . " 
else 

rm  -f  STATS. prt 

touch  STATS. prt  1 

echo  "Brawler  Post-Processor"  >>  STATS. prt 

echo  "  Execution  Date /Time:  $CurrentDate  /  $CurrentTime"  »  STATS. prt 
echo  "  "  »  STATS. prt 

echo  "  Task:  $Task"  »  STATS. prt 

echo  "  Directory:  'basename  $cwd' "  »  STATS. prt 

echo  "  "  >>  STATS. prt 

# 

echo  "  Task:  $Task"  »  $Terminal 

echo  "  Directory:  'basename  $cwd' "  >>  $Terminal 

echo  "  "  »  $Terminal 

# 

set  InputFiles=" ' Is  LOG.*  |  sort  -t.  +ln' " 

# 

Is  database  >&  /dev/null 
set  DBs tatus=$ status 

if  ( $DBstatus  ==  2  |  $OverWrite  ==  $Yes)  then 
echo  "  Creating  Database  File. . . " 

rm  -f  database 
touch  database 
# 

if  (-e  STATS.dat)  then 

nawk  -f  $RunPath/read__input . awk  STATS.dat  >>  database 
endif 
# 

foreach  InputFile  ($InputFiles) 

echo  -n  ’  Processing  LOG  File  " • $InputFile’ " . . . ' 

rm  -f  database. tmp 

nawk  -v  PrePro=$PreProcess  -f  $RunPath/build_database . awk  $InputFile 
>  database . tmp 

if  ($status  ==  1)  then 
echo  "  (Good) " 
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cat  database. tmp  »  database 
else 

echo  H  (Bad) " 

set  Field2=' echo  $InputFile  |  cut  -f2  -d.' 
mv  $InputFile  log.$Field2 
end  if 

end 

rm  -f  database . tmp 
endif 
# 

#  echo  "  Generating  Aircraft  Survivability  Report..." 

#  nawk  -f  $RunPath/ac_summary . awk  database  >>  STATS. prt 

# 

##  echo  "  Generating  Aircraft  Kill  Report..." 

##  /usr/ local /bin /perl  $RunPath/kill .perl  database  LOG.l  >>  STATS. prt 

# 

#  echo  "  Generating  Mission  Success  Report..." 

#  nawk  -f  $RunPath/msn_success .awk  database  >>  STATS. prt 

.  ft _ _ _ 

echo  "  Generating  Missile  Effectiveness  Report..." 

nawk  -f  $RunPath/msl_summary.awk  database  »  STATS. prt 

_ # _ _ _ 

#  echo  "  Generating  Missile  Launch  Report..." 

#  nawk  -f  $RunPath/msl_launch . awk  database  »  STATS. prt 

# 

#  echo  "  Generating  Missile  Launch  Distribution  Reports..." 

#  nawk  -f  $RunPath/msl_launch_dis_bvr . awk  database  >>  STATS. prt 

#  nawk  -f  $RunPath/msl_launch_dis_wvr .awk  database  »  STATS. prt 

# 

#  echo  "  Generating  Missile  Shot/Kill  Breakdown  Report..." 

#  nawk  -f  $RunPath/msl_breakdown.awk  database  »  STATS. prt 

# 

#  echo  "  Generating  Missile  Range/Aspect  Report..." 

#  nawk  -f  $RunPath/msl_rng_asp . awk  database  >>  STATS. prt 

# 

#  echo  "  Generating  Gun  Reports ..." 

#  nawk  -f  $RunPath/gun_rng . awk  database  »  STATS. prt 

#  nawk  “f  $RunPath/gun_aspect .awk  database  »  STATS. prt 

#  nawk  -f  $RunPath/gun_rng_aspect . awk  database  »  STATS. prt 

#  nawk  -f  $RunPath/gun_rng_aspect_pk . awk  database  »  STATS. prt 

#  nawk  -f  $RunPath/gun_time_rng . awk  database  >>  STATS. prt 

# 

#  echo  "  Generating  Weapon  Time/Range  Report..." 

#  nawk  -f  $RunPath/wpn_time_rng . awk  database  »  STATS. prt 
# 

endif 

end 

# 

exit 
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Appendix  D.  STATS.ort  Missile  Effectiveness  Reports 


MISSILE  EFFECTIVENESS  REPORT 


MSLI 

MSLR 

Firings 

537 

1759 

Fuzings 

294 

1108 

Kills 

172 

626 

Kills/Firing 

32% 

35% 

BLUE: 

RF  Msl  Firings/Kills 

1759 

626 

IR  Msl  Firings /Kills 

483 

165 

RED: 

RF  Msl  Firings/Kills 

0 

0 

IR  Msl  Firings/Kills 

54 

7 

THUNDER  MISSILE  EFFECTIVENESS  REPORT 

MSLI  MSLR 

Adj .  Kills/Firing  0.32  0.36 
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Author:  LAWRENCE  J.  TARANTO 

ASC/XRA,  BUILDING  16 
2275  D  STREET  SUITE  10 
WRIGHT- PATTERSON  AFB  OH  45433-7227 

Purpose:  To  generate  missile  summary  report. 

Date:  15  OCTOBER  1997 

BEGIN  { 

FS="  |  ■ 

cases=0 
b_total=0 
r_total=0  } 

# 

{ 

rec__type=$2 

if  (rec_type  ==  "MSL" )  { 

mslid=$3 
ms typ=$4 

launcher [mstyp] =$5 
freason=$10 
fuzee=$ll 
kill=$12 

if  (kill  ==  "T")  {kill=l} 
if  (kill  ==  "F")  {kill=0} 

# 

if  (f reason  !=  "")  { 

failure [f reason, mstyp] ++ 
fmode [ treason] ++ 
failedfmstyp] ++ 

if  (fuzee  !=  0)  {  fuz_n__f ail  [  treason, mstyp] ++  } 

} 

# 

mode=$14 

launch_mode  [mstyp, mode]  ++ 

if  (kill  ==  1)  {  lethal_launch_mode [mstyp , mode] ++  } 

b_init=$15 

shooter=int ( substr (mslid,  1,  2)) 
if  (shooter  >  b_init)  {  blue_red=2  } 
if  (shooter  <=  b_init)  {  blue_red=l  } 
if  (mode  ==  » SEMI„ACTIVE_SEEKER_LOCK" )  { 

modes [mode] ++ 
rf_f ires [blue_red] ++ 

if  (kill  ==  1)  {  rf_kills [blue_red] ++  }  } 

else  if  (mode  ==  "CMD_GUIDED_DES_RDR" )  { 

modes [mode] ++ 
rf_f ires [blue_red] ++ 

if  (kill  ==  1)  {  rf_kills [blue_red] ++  }  } 

else  if  (mode  ==  " CMD_GUIDED_DES_ITB " )  { 

modes [mode] ++ 
rf_f ires [blue_red] ++ 

if  (kill  ==  1)  {  rf_kills [blue_red] ++  }  } 

else  if  (mode  ==  " CMD__GUIDED_DES_OTD " )  { 

modes [mode] ++ 
rf_f ires [blue_red] ++ 

if  (kill  ==  1)  {  rf__kills  [blue_red]  ++  }  } 
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else  if  (mode  ==  "CMD_GUIDED_DES„STT" )  { 

modes [mode] ++ 
rf_f ires [blue_red] ++ 

if  (kill  ==  1)  {  rf_Jcills  [blue_red]  ++  }  } 
else  if  (mode  ==  "CMD_GUIDED_UNDES_VIS " )  { 
modes [mode] ++ 
rf_f ires [blue__red] ++ 

if  (kill  ==  1)  {  rf_kills [blue_red] ++  }  } 

else  if  (mode  ==  "CMD_GUIDED_UNDES_IRST" )  { 
modes [mode] ++ 
rf_f ires [blue_red] ++ 

if  (kill  ==  1)  {  rf_kills [blue„red] ++  }  } 
else  if  (mode  ==  " PASSIVE_SEEKER_LOCK_IR" )  { 
modes [mode] ++ 
ir_f ires [blue_red] ++ 

if  (kill  ==  1)  {  ir_kills [blue_red] ++  }  } 
else  if  (mode  ==  " LOCK_ON_AFTER_LN_HMS ” )  { 

modes [mode] ++ 
ir_f ires [blue_red] ++ 

if  (kill  —  1)  {  ir_kills [blue_red] ++  }  } 
else  if  (mode  ==  "LOCK_ON_AFTER_LN_RDR" )  { 

modes [mode] ++ 
ir_f ires [blue_red] ++ 

if  (kill  ==  1)  {  ir_kills [blue_red] ++  }  } 
else  { 

printf  (  "  \n%s\n"  ,  "»>  WARNING:  UNKNOWN  MISSILE  LAUNCH  MODE 
ENCOUNTERED  «<")  } 

# 

launchings [mstyp] ++ 

if  (freason  ==  "  FUZED  ON  DEAD  TARGET  ")  {  fuzed^dead [mstyp] ++  } 

if  (fuzee  !=  0)  {  fuzed [mstyp] ++  } 

if  (kill  ==  1)  {  killed [mstyp] ++  } 

} 

if  (rec_type  ==  "AC  " )  { 

cases=cases+$3 
b_total=b_total+$5 
r_total=r_total+$ll 

} 

} 

# 

END  { 

# _ _ 

#printf  " \n\n%27s " ,  "MISSILE  EFFECTIVENESS  REPORT\n" 

#printf  " % 1 0 0 s " , 

=====================“=:======\n\nn 

# 

#printf ( " %25s " ,  "  " ) 

#for  (i  in  launchings)  {  printf (" %13s ", i)  } 

#printf ( " \n%25s " , "  " ) 

#for  (i  in  launchings)  {  print f  ( "  - " )  } 

#printf ( " \n%-19s" , "Firings  " ) 

#f or  (i  in  launchings)  {  printf ( " %13d" , launchings [i] )  } 

#printf ( " \n%-19s " , " Fuzings  " ) 

#for  (i  in  launchings) {  printf (" %13d" f fuzed [ i ] )  } 

#printf ( " \n%-19s" , "Kills  " ) 

#for  (i  in  launchings)  {  printf ( n%13d" , killed [i] )  } 

#printf ( " \n%-25s " , "Fuzings/Firing" ) 
for  (i  in  launchings)  { 

rate=100 . 0* (fuzed[i] /launchings [i] ) 

#  printf ("%7d%%%5s'\ rate,  "  " ) 
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} 

#printf ( " \n%-25s" , "Kills/Firing  " ) 
for  (i  in  launchings)  { 

rate=100 . 0* (killed[i] /launchings [i] ) 

#  printf (" %7d%%%5s ",  rate,  "  " ) 

} 

#printf ( ” \n%-25s" , "Kills/Fuzing  " ) 
for  (i  in  launchings)  { 
rate=0 . 0 

if  ( fuzed [ i ] !=0)  rate=100.0* (killed [i] / fuzed [i] ) 

#  printf ( " %7d%%%5s " , rate , "  " ) 

} 

tprintf ( " \n" ) 

#printf ( " \n%-5s" , "BLUE: " ) 

#printf ( " \n%-25s" , n  RF  Msl  Firings/Kills”) 

#printf  ( ,,%13d%13d" ,  rf„fires  [1]  ,  rf_kills  [1]  ) 

#printf ( " \n%-25s" , "  IR  Msl  Firings/Kills " ) 

#printf ( " %13d%13d" , ir_fires [1] , ir_kills [1] ) 

#printf ( " \n%-4s " , "RED: " ) 

#printf ( " \n%-25s " , "  RF  Msl  Firings/Kills") 

#printf ( " %13d%13d" , rf_f ires [2] , rf_kills [2 ] ) 

#printf ( " \n%-25s " , "  IR  Msl  Firings/Kills " ) 
tprintf ( " %13d%13d" , ir_fires [2] , ir_ kills [2] ) 

# 

#printf  "\n\n\n%36s",  "THUNDER  MISSILE  EFFECTIVENESS  REPORT \n" 

#printf  "%100s", 

=  =  =  =  =  =  =  =  =  =  “  =  =  :=  =  =  =  =  =  =  =  =  =  =  =  =  =  = \n\n" 

# 

#printf ( " %25s " , "  " ) 

#for  (i  in  launchings)  {  printf (" %13s ", i )  } 

#printf ( " \n%25s " , "  ") 

#for  (i  in  launchings)  {  printf  ( "  - ")  } 

printf ( " \n%-25s " , "Adj .  Kills/Firing  " ) 
for  (i  in  launchings)  { 
if  ( fuzed [i]  >  0)  { 

den= launchings [ i ] - ( fuzed_dead [ i ] * ( launchings [ i ] / fuzed [ i ] ) ) 
rate=killed [ i ] /den  } 
else  {rate-0.0} 

printf ( " %7 . 2f %6s " , rate, "  " ) 

} 

# 

printf ( " \n\n" ) 
b_f ires=0 
r_f ires=0 

for  (i  in  launchings)  { 

if  (launcher [i]  <=  b_total /cases )  {  b_fires=b_f ires+launchings [ i ]  } 
if  (launcher [i]  >  b_total/cases )  {  r_fires=r_f ires+launchings [ i ]  } 

} 

tprintf  " % 2 5 s  %6 . 2 f \n" , "BLUE  Shots/Aircraft : " ,  b_f ires/b_total 
#printf  "%25s  %6.2f","RED  Shots/Aircraft :" ,  r_f ires/r_total 

# _ _ _ 

printf  " \n\n\n%23s" ,  "MISSILE  FAILURES  RE PORT \n"  ~~ 

printf  "%100s". 


#  XnXn 
printf ( " %-25s " , "  FAILURE  MODE" ) 
for  (i  in  launchings)  {  printf (" %13s ", i)  } 
printf ( " \n%25s" , "  " ) 
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for  (i  in  launchings)  {  printf ( "  - " )  } 

printf ( " \n" ) 
for  (i  in  fmode)  { 
printf { "%-25s" , i) 
for  (1  in  launchings)  { 
row[l, 1] = failure [i , 1] 
row[l, 2] =fuz_n_fail [1,1] 
sum[l, 1] +=failure[i, 1] 
sum [ 1 , 2] +=fuz_n_fail [i, 1]  } 

for ( 1  in  launchings)  {  printf ( " %7d/%5d” , row [ 1 , 1 ], row [ 1 , 2 ] )  } 

printf ( " \n" )  } 

# 

printf ( " \n%-25s" , "  FAILED  BUT  FUZED") 

for  (1  in  launchings)  {  printf ( " %7d%6s ", sum[ 1 , 2 ], "  ")  } 
printf ( "\n%-25s" , "  HARD  FAILURES") 

for  (1  in  launchings)  {  printf ( " %7d%6s " , sum [ 1 , 1 ] -sum[l , 2 ] , "  ")  } 

# 

printf  " \n\n\n%26s " ,  "MISSILE  LAUNCH  MODE  REPORT \n" 
printf  "%100s". 


======================:=======  \.n\n" 

# 

printf ( " %-25s " , "  LAUNCH  MODE") 

for  (i  in  launchings)  {  printf (" %13s ", i )  } 

printf ( " \n%25s " , "  " ) 

for  (i  in  launchings)  {  printf  ( "  - " )  } 

printf ( " \n" ) 
for  (i  in  modes)  { 
printf ( "%-25s" , i) 

for  (1  in  launchings)  {  printf (" %7d% 6s ", launch_mode [1 , i] , "  ")  } 

printf ("\n")  } 

# 

printf  " \n\n\n%26s " ,  "MISSILE  LETHAL  LAUNCH  MODE  REPORT\n" 
printf  " %100s " , 


======================== =====\n\n" 

# 

printf ( "%-25s" , "  LETHAL  LAUNCH  MODE") 
for  (i  in  launchings)  {  printf (" %13s ", i)  } 
printf ( " \n%25s" , "  " ) 

for  (i  in  launchings)  {  print f  ( "  - ")  } 

printf ("\n") 
for  (i  in  modes)  { 
printf ("%-25s" ,i) 

for  (1  in  launchings)  {  printf {" %7d% 6s ", lethal_launch_mode [ 1 , i ] 

} 

printf ("\n")  }  } 

} 
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Appendix  F.  make  flyAKPF 


#  C- she 11  program  make_flyAKPF 

# 

#  Author:  Bryan  Livergood 

# 

#  This  program  searches  all  BRAWLER  log  files  in  a  directory 

#  and  obtains  Firings  and  Adjusted  Kills  per  Firing  (AKPF) 

#  data  for  each  log.  It  then  creates  one  file  called 
fly.FAKPF 

#  that  contains  a  list  of  all  the  Firings  values  followed  by 
a 

#list  of  all  the  AKPF  values  for  each  BRAWLER  run. 

# 

#!/bin/csh 
mkdir  akpfSTATS 
foreach  fyle  ('Is  LOG.*  $1') 
echo  'processing  '$fyle 
grep  ALIVE  $fyle 
echo  'grep  status  ' $status 
if  (!($status))  then 

set  fylenum  =  'Is  $fyle  |  cut  -dG  -f2' 

echo  'log  number  '$ fylenum 

mkdir  temporary 

cp  $fyle  . /temporary 

cd  temporary 

STATS 

grep  -h  'Firings'  STATS. prt  »  . . / f ly . Fir e+AKPF 
grep  -h  'Adj.  Kills/Firing'  STATS. prt  » 

. . / fly . Fire+AKPF 

mv  STATS. prt  . . /akpf STATS/STATS$fylenum 
rm  * 
cd  .  . 

rmdir  temporary 
endif 

end 

grep  -h  'Adj.  Kills/Firing'  fly . Fire+AKPF  >  fly.FAKPF 
grep  -h  'Firings  ’  fly .Fire+AKPF  »  fly.FAKPF 
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Appendix  G.  BRAWLER  AKPF  and  Firings  Data 


Replication 

Number 


1 


Scenario  1 

Scenario  2 

Scenario  3 

AKPF 

Firings 

AKPF 

Firings 

AKPF 

Firings 

0.40 

5 

0.00 

6 

0.0 

0 

0.00 

4 

0.00 

0 

0.0 

2 

0.25 

4 

0.00 

5 

0.4 

5 

0.00 

1 

0.00 

0 

0.5 

2 

0.50 

4 

0.67 

3 

0.7 

3 

0.20 

5 

0.12 

8 

1.0 

1 

0.25 

4 

0.33 

3 

0.0 

3 

0.50 

4 

0.50 

2 

0.0 

1 

0.00 

2 

0.67 

3 

0.2 

6 

1.00 

1 

0.25 

4 

0.7 

3 

0.00 

0 

0.50 

2 

0.0 

0 

0.00 

2 

0.33 

3 

0.5 

2 

0.00 

1 

0.00 

1 

0.0 

0 

0.20 

5 

0.00 

2 

0.3 

3 

0.00 

0 

0.40 

5 

0.5 

2 

0.00 

2 

0.50 

2 

0.0 

4 

0.50 

2 

1.00 

2 

0.3 

4 

0.40 

5 

0.00 

4 

0.0 

2 

0.00 

0 

0.00 

2 

0.0 

1 

0.33 

3 

0.00 

3 

1.0 

2 

1.00 

2 

1.00 

1 

1.0 

1 

1.00 

1 

0.20 

5 

0.5 

2 

0.00 

1 

0.00 

2 

0.0 

0 

0.00 

0 

0.25 

4 

0.0 

0 

0.67 

3 

0.67 

3 

0.5 

2 

1.00 

2 

0.00 

3 

0.7 

3 

0.00 

5 

0.33 

3 

0.0 

2 

1.00 

1 

0.67 

3 

0.0 

1 

0.50 

2 

0.00 

2 

0.0 

0 

0.00 

2 

0.00 

5 

0.0 

1 

0.00 

0 

0.67 

3 

1.0 

2 

0.00 

3 

0.25 

4 

0.0 

0 

0.00 

0 

0.33 

3 

0.0 

2 

0.25 

4 

1.00 

1 

0.5 

2 

0.00 

0 

0.33 

3 

0.0 

2 

0.00 

3 

0.33 

3 

0.0 

0 

0.00 

4 

0.50 

2 

1.0 

1 

1.00 

1 

0.00 

1 

1.0 

2 

0.33 

3 

1.00 

1 

0.0 

0 

1.00 

1 

0.17 

6 

1.0 

1 

Replication 

Scenario  1 

j  Scenario  2 

Scenario  3 

Number 

AKPF 

Firings 

AKPF 

Firings 

AKPF 

Firings 

41 

0.33 

3 

■B| 

0 

4 

42 

0.20 

5 

ESI 

0 

0.0 

4 

43 

0.00 

3 

ESI 

1 

0.5 

2 

44 

1.00 

1 

ESI 

4 

0.0 

4 

45 

0.00 

1 

0.00 

1 

0.5 

2 

46 

0.17 

6 

1.00 

1 

0.0 

2 

47 

0.00 

1 

0.25 

4 

0.3 

3 

48 

0.00 

0 

0.50 

2 

1.0 

1 

49 

0.17 

6 

1.00 

1 

0.0 

0 

50 

0.50 

2 

0.50 

2 

0.3 

6 

51 

0.17 

6 

0.00 

0 

1.0 

1 

52 

0.00 

0 

0.33 

3 

0.0 

1 

53 

0.20 

5 

0.17 

6 

0.3 

3 

54 

0.00 

0 

0.00 

0 

0.0 

2 

55 

0.33 

3 

0.00 

1 

0.0 

0 

56 

0.50 

4 

0.33 

3 

0.7 

3 

57 

0.00 

1 

0.50 

2 

0.5 

2 

58 

0.50 

4 

0.50 

2 

0.0 

8 

59 

1.00 

2 

1.00 

1 

0.0 

0 

60 

0.00 

3 

0.00 

2 

0.0 

1 

61 

0.00 

2 

0.00 

3 

1.0 

1 

62 

0.00 

0 

0.50 

4 

0.0 

1 

63 

0.00 

0 

0.00 

0 

0.3 

3 

64 

0.50 

2 

0.00 

3 

0.0 

4 

65 

0.00 

0 

0.50 

4 

0.0 

2 

66 

1.00 

2 

0.00 

2 

1.0 

2 

67 

0.00 

0 

0.00 

2 

0.1 

7 

68 

1.00 

2 

0.00 

3 

0.0 

4 

69 

0.00 

3 

0.00 

3 

0.0 

1 

70 

0.00 

4 

0.40 

5 

0.0 

0 

71 

0.50 

2 

0.00 

0 

0.0 

1 

72 

0.20 

5 

1.00 

1 

0.0 

0 

73 

0.00 

2 

0.20 

5 

0.0 

1 

74 

1.00 

1 

0.00 

3 

0.0 

0 

75 

0.40 

5 

0.25 

4 

0.0 

2 

76 

0.00 

2 

0.00 

0 

0.0 

1 

77 

0.25 

4 

0.00 

1 

1.0 

1 

78 

1.00 

1 

0.50 

2 

0.0 

2 

79 

1.00 

1 

0.00 

5 

0.0 

1 

80 

0.00 

0 

0.00 

1 

0.0 

1 
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Replication 

Scenario  1 

Scenario  2 

|  Scenario  3 

Number 

AKPF 

Firings 

AKPF 

Firings 

AKPF 

Firings 

81 

0.00 

0 

0.00 

1 

LO 

o 

2 

82 

0.33 

3 

0.43 

7 

0.0 

1 

83 

1.00 

1 

0.40 

5 

1.0 

1 

84 

0.00 

1 

0.50 

2 

0.5 

2 

85 

0.00 

2 

0.00 

3 

0.4 

5 

86 

0.33 

3 

1.00 

3 

0.0 

1 

87 

0.50 

2 

0.00 

2 

0.0 

0 

88 

0.50 

4 

0.50 

2 

0.5 

2 

89 

1.00 

1 

0.00 

3 

0.0 

3 

90 

0.00 

0 

0.00 

0 

0.2 

6 

91 

1.00 

2 

0.00 

0 

0.0 

1 

92 

1.00 

1 

0.00 

1 

0.0 

1 

93 

0.00 

0 

0.33 

3 

0.3 

4 

94 

1.00 

1 

0.17 

6 

0.0 

1 

95 

0.17 

6 

0.00 

2 

0.3 

3 

96 

1.00 

2 

0.50 

2 

0.0 

0 

97 

0.14 

7 

0.33 

3 

0.0 

0 

98 

0.50 

2 

0.20 

5 

0.5 

4 

99 

0.00 

0 

0.00 

2 

0.3 

3 

100 

0.67 

3 

0.33 

3 

0.0 

1 

101 

0.00 

0 

0.50 

2 

0.0 

0 

102 

1.00 

1 

1.00 

1 

0.0 

1 

103 

0.00 

4 

0.33 

3 

0.0 

4 

104 

0.00 

2 

0.00 

1 

0.0 

2 

105 

0.33 

3 

0.25 

4 

1.0 

1 

106 

0.14 

7 

0.33 

3 

0.3 

3 

107 

0.00 

5 

0.00 

5 

0.0 

6 

108 

0.00 

0 

0.50 

2 

0.0 

0 

109 

0.00 

2 

1.00 

2 

1.0 

1 

110 

0.00 

0 

0.00 

1 

1.0 

1 

111 

0.50 

2 

0.00 

0 

1.0 

1 

112 

0.00 

2 

0.00 

4 

0.5 

4 

113 

0.00 

0 

0.40 

5 

0.3 

3 

114 

0.50 

4 

1.00 

1 

0.0 

0 

115 

0.00 

1 

0.25 

4 

0.4 

5 

116 

0.00 

0 

0.33 

3 

0.5 

2 

117 

0.00 

3 

0.50 

4 

0.0 

1 

118 

0.50 

2 

0.25 

4 

0.0 

3 

119 

0.00 

0 

0.20 

5 

0.3 

3 

120 

0.00 

0 

0.00 

4 

1.0 

1 
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Replication 

Scenario  1 

Scenario  2 

Scenario  3 

Number 

AKPF 

Firings 

AKPF 

Firings 

AKPF 

Firings 

mSm 

0.33 

3 

1.00 

1 

0.0 

1 

apl 

0.50 

4 

0.00 

2 

0.1 

7 

0.10 

10 

0.50 

2 

0.2 

5 

124 

0.33 

3 

0.25 

4 

0.0 

0 

125 

0.00 

0 

1.00 

1 

0.0 

0 

126 

0.00 

0 

0.00 

2 

0.1 

9 

127 

0.00 

4 

0.50 

2 

0.0 

0 

128 

0.00 

7 

0.67 

3 

0.3 

4 

129 

0.33 

3 

0.33 

3 

1.0 

i 

130 

0.00 

0 

0.00 

2 

0.0 

0 

131 

0.33 

6 

0.00 

1 

0.0 

1 

132 

1.00 

1 

0.33 

3 

0.0 

1 

133 

0.00 

2 

0.25 

4 

0.0 

0 

134 

0.00 

1 

0.00 

1 

1.0 

1 

135 

0.50 

2 

0.40 

5 

0.3 

3 

136 

0.25 

4 

0.60 

5 

0.0 

2 

137 

0.20 

5 

0.33 

3 

0.0 

2 

138 

1.00 

3 

0.67 

3 

0.0 

0 

139 

0.67 

3 

1.00 

1 

0.0 

0 

140 

0.00 

1 

0.00 

1 

0.0 

0 

141 

0.00 

1 

0.00 

1 

0.0 

0 

142 

0.00 

0 

0.25 

4 

0.0 

0 

143 

0.00 

3 

0.00 

0 

0.0 

2 

144 

0.50 

2 

0.33 

3 

0.0 

0 

145 

0.00 

0 

0.33 

3 

0.0 

4 

146 

0.25 

4 

1.00 

3 

1.0 

1 

147 

1.00 

1 

0.33 

3 

0.0 

0 

148 

0.50 

2 

0.50 

4 

0.3 

3 

149 

0.50 

2 

0.00 

2 

0.4 

5 

150 

0.33 

3 

0.00 

0 

0.3 

4 

151 

0.00 

0 

0.33 

3 

0.0 

0 

152 

0.00 

0 

0.00 

0 

0.0 

0 

153 

0.50 

2 

1.00 

1 

0.3 

3 

154 

0.25 

8 

0.50 

2 

0.3 

6 

155 

1.00 

1 

1.00 

1 

0.0 

0 

156 

0.33 

3 

1.00 

1 

0.5 

2 

157 

0.00 

0 

0.50 

2 

0.0 

0 

158 

0.00 

1 

0.50 

2 

0.0 

0 

159 

0.00 

1 

0.67 

3 

1.0 

2 

160 

0.75 

4 

0.20 

5 

0.0 

3 
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Replication 

Scenario  1 

Scenario  2 

|  Scenario  3 

Number 

AKPF 

Firings 

AKPF 

Firings 

AKPF 

Firings 

mm 

0.00 

3 

1 

0.0 

3 

0.00 

0 

0.5 

4 

Hgfl 

0.00 

0 

1.00 

1 

0.0 

3 

164 

0.00 

1 

1.00 

1 

0.0 

0 

165 

0.00 

3 

0.00 

1 

0.0 

0 

166 

1.00 

1 

0.00 

3 

0.0 

1 

167 

0.50 

2 

0.40 

5 

0.3 

3 

168 

0.00 

1 

0.50 

2 

0.0 

1 

169 

0.50 

2 

0.00 

0 

0.4 

5 

170 

1.00 

1 

1.00 

1 

0.0 

2 

171 

0.50 

4 

0.00 

4 

1.0 

1 

172 

0.00 

1 

0.00 

0 

0.0 

2 

173 

0.50 

2 

0.50 

2 

0.0 

3 

174 

0.00 

2 

0.67 

3 

0.0 

0 

175 

0.00 

0 

0.25 

4 

1.0 

1 

176 

1.00 

1 

0.50 

4 

0.0 

0 

177 

0.40 

5 

0.50 

2 

0.0 

0 

178 

0.43 

7 

0.00 

0 

0.0 

1 

179 

0.00 

3 

0.33 

3 

0.5 

4 

180 

0.00 

2 

0.50 

2 

1.0 

1 

181 

0.14 

7 

0.00 

4 

0.5 

2 

182 

0.50 

6 

0.00 

0 

0.3 

4 

183 

0.00 

0 

0.50 

6 

0.0 

2 

184 
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Appendix  H.  Derivation  of  Adjusted  Kills  per  Firing  (AKPF) 


The  following  equations  are  presented  in  Section  4.4.1  and  are  presented  again 

here: 


KPF  = 


E Kills 
Firings 


(4.2) 


AKPF  = 


E  Kills 


y  Firin  S  Fuzin8sDEAD  '  X  Firin8s 
w  6  E Fuzings 


(4.3) 


AKPF  =  KPF  ■ 


E Fuzings 

£  Fuzings  -  E  Fuzings  DEAD ) 


(4.4) 


These  are  the  equations  that  STATS  uses  to  calculate  kills  per  firing  (KPF)  and 
adjusted  kills  per  firing  (AKPF). 

We  now  present  the  derivation  of  (4.4)  from  (4.3).  The  actual  computation 
performed  by  STATS  follows: 


AKPF  = 


^  Kills 


where 


den  =  E Firings 


^ Fuzings DEAD  -  ^Firings 
E Fuzings 


H-l 


which  we  manipulate  as  follows: 


den  _  X  FirinSs '  X  Fuzings  X  FuzinSsDEAD  '  X  Firings 
X Fuzings  ^ Fuzings 


to  obtain: 


^  _  ^Firings  •  Fuzings  -  ^Fuzings DEAD ) 

X Fuzings 


Since 


^  Kills 


AKPF  =  -  =>  AKPF  =  Y  Kills  •  — 

den  den 


This  implies: 


AKPF  =  \  Kills  ■  -= -  Yl  Fuzmgs  - 

2j  Firin8S  •  12,  Fuzings  -  2,  Fuzings DEAD  j 

X  Kills  Fuzings 

Yj  Firings-  (X  Fuzings  -  Fuzings  DEAD  ) 


which  gives  the  desired  result: 


AKPF  =  KPF  • 


^ Fuzings 

(X  Fuzings  -  £  Fuzings  DEAD ) 


(4.4) 
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